1 Introduction↩︎

Binarized neural networks (BNNs) are feedforward neural networks with binary weights and activation functions ([1]). BNNs consist of the input layer \(\ell=0\), \(L\) hidden layers \(\ell=1,\ldots,L\), and the output layer \(\ell=L+1\). Each layer \(\ell\in \{0,1,\ldots,L+1\}\) consists of \(n^\ell\) neurons, \(N^\ell = \{1, \ldots, n^\ell \}\). Every \(i\in N^\ell\) is connected to every \(j\in N^{\ell- 1}\) with weight \(W_{ij}^\ell\in \{-1, 0, 1\}\) and has bias \(b_i^\ell\in \mathbb{Z}\) for \(\ell\in \{1, \ldots, L+ 1\}\).

Thanks to binary weights and activation functions, BNNs drastically reduce memory size and accesses and substantially improve power-efficiency over standard neural networks, replacing most arithmetic operations with bit-wise operations ([1]). BNNs have also achieved nearly state-of-the-art results in image classification ([1]). Moreover, BNNs have achieved results comparable to feedforward neural networks in image detection ([2]), image super resolution ([3]), and text classification ([4]). For these reasons, BNNs have been applied in small embedded devices ([5]).

In this work, we study the verification problem associated with a given BNN that is used to classify feature vectors. BNNs map a feature vector in \([0, 1]^{n^0}\) whose coordinates are quantized as multiples of \(\frac{1}{q}\) for \(q\in \mathbb{N}\) to the real output vector by the function \(f= (f_1, \ldots, f_{n^{L+ 1}}): \frac{1}{q} \mathbb{Z}_+^{n^0} \cap [0, 1]^{n^0} \to \mathbb{R}^{n^{L+ 1}}\) which is defined recursively using weights and biases in each layer. (See Section 2 for the detailed definition of \(f\).) Each \(t\in N^{L+ 1}\) corresponds to a class. A feature vector \(\bar{\mathbf{x}}\in \frac{1}{q} \mathbb{Z}_+^{n^0} \cap [0, 1]^{n^0}\) is classified as the class \(\bar{t}\in N^{L+ 1}\) that corresponds to \(\mathop{\mathrm{arg\,max}}_{t\in N^{L+ 1}} f_t(\bar{\mathbf{x}})\).

Verifying BNNs against input perturbation for multiple feature vectors provides a measure of their robustness. The BNN verification problem is defined for a given trained BNN, a given feature vector \(\bar{\mathbf{x}}\) with class \(\bar{t}\), and input perturbation \(\epsilon> 0\). We say that \(\bar{\mathbf{x}}\) is \(\epsilon\)-verified if the BNN classifies \(\mathbf{x}^0\) as \(\bar{t}\) for every feature vector \(\mathbf{x}^0 \in \frac{1}{q} \mathbb{Z}_+^{n^0} \cap [0, 1]^{n^0}\) having distance between \(\mathbf{x}^0\) and \(\bar{\mathbf{x}}\) at most \(\epsilon\). The BNN verification problem is to determine whether a given \(\bar{\mathbf{x}}\) can be \(\epsilon\)-verified in the BNN. The answer is true if and only if there does not exist a feature vector \(\mathbf{x}^0\) with distance at most \(\epsilon\) from \(\bar{\mathbf{x}}\) satisfying \(f_t(\mathbf{x}^0) > f_{\bar{t}}(\mathbf{x}^0)\) for some class \(t\ne \bar{t}\), which is equivalent to the following maximum \(z^*_\epsilon(\bar{\mathbf{x}})\) being non-positive: \[z^*_\epsilon(\bar{\mathbf{x}}) = \max_{\substack{\mathbf{x}^0 \in \frac{1}{q} \mathbb{Z}_+^{n^0} \cap [0, 1]^{n^0}\\ \lVert\mathbf{x}^0 - \bar{\mathbf{x}}\rVert_p \le \epsilon}} \Big\{f_t(\mathbf{x}^0) - f_{\bar{t}}(\mathbf{x}^0): t\in N^{L+ 1} \setminus \{\bar{t}\}\Big\}.\label{eq:bnnVeri}\tag{1}\] Here, the maximum perturbation from \(\bar{\mathbf{x}}\) is defined using an \(\ell_p\) norm. For the most part, in this work we focus on the case of \(p= 1\), but the integer programming (IP) approaches we consider can also be directly applied to the cases of \(p= \infty\) and \(p= 2\), with the caveat that when using an \(\ell_2\) norm, the resulting IP formulations we consider become integer (convex) quadratic programs, as opposed to integer linear programs.

This paper contributes to the literature exploring IP methods to solve the BNN verification problem. IP methods to solve the BNN verification problem can be implemented by solving 1 to optimality – we refer to this as the BNN verification optimization problem. However, to verify a BNN, the process of solving 1 can be terminated as soon as the sign of \(z^*_\epsilon(\bar{\mathbf{x}})\) is determined. Specifically, the solution process of 1 can be terminated when a feasible solution with a positive objective value is found, or a non-positive upper bound for the optimal objective value is obtained. If we terminate an optimization algorithm for solving 1 according to this condition, we refer to this as the BNN verification problem.

Our first contribution is to present a new way for obtaining a linear objective of the BNN verification problem by creating a single optimization problem that incorporates the decision of which alternative class a perturbed feature vector is misclassified as, rather than considering each alternative class separately. Our second contribution is to describe a new technique for generating valid inequalities for the IP formulation by exploiting the recursive structure of BNNs. These valid inequalities, called layerwise derived valid inequalities, are generated by considering one layer at a time and solving IP subproblems to check validity among a set of natural candidate inequalities. While deriving these valid inequalities requires solving IP subproblems, these IP subproblems are much easier to solve than the full BNN verification problem because they consider only a single layer at a time and do not include any of the “big-\(M\)” constraints required for the IP formulation of the BNN verification problem.

[6] study the problem of attacking a BNN, which is equivalent to the BNN verification problem. They present a mixed-integer linear programming (MILP) formulation of the problem and propose a heuristic for generating solutions that works by propagating a target from the last layer to the first layer by solving a single-layer MILP at each step. This approach is proven to be successful at finding successful attacks (i.e., showing an input \(\bar{x}\) cannot be \(\epsilon\)-verified) but cannot positively verify a solution. Interestingly, our method also proceeds in layers, but from the initial layer forwards, and its strength is in being able to \(\epsilon\)-verify inputs. [7] also study the MILP formulation of the BNN verification problem and study the convex hull for a single neuron to improve this MILP formulation. [8] extended Han and Gómez’s work by exploring the convex hull for a pair of neurons in the same hidden layer. Using this convex hull, they proposed a constraint generation framework to solve the MILP problem for the BNN verification problem. In related work on feedforward (not binarized) neural networks, [9] addressed the problem of verifying such networks. They formulated this problem as an MILP problem and exploited the idea of fixing variables to solve this MILP problem. Our layerwise derived valid inequalities represent an adaptation and extension of this approach to the BNN setting.

Another approach for solving the BNN verification problem is to formulate it as a Boolean satisfiability problem ([10]–[13]). These papers formulated the condition that \(\bar{\mathbf{x}}\) can be \(\epsilon\)-verified in a given BNN with Boolean formulas in the case that the maximum perturbation from \(\bar{\mathbf{x}}\) is defined using an \(\ell_\infty\) norm. [14] explored a different method for solving the BNN verification problem by describing the set of all feasible output vectors or an overapproximation of this set as a star set and checking whether this star set contains an output vector from a perturbed feature vector misclassified by the BNN. This work also handles the case that the maximum perturbation is defined using an \(\ell_\infty\) norm. [15] and [16] addressed a different problem of counting the number of perturbed feature vectors misclassified by the BNN in the case that the maximum perturbation is defined using an \(\ell_1\) norm. We focus on the IP approach due to its inherent flexibility, for example, in easily adapting the set of allowed perturbations from \(\bar{\mathbf{x}}\) to be defined by an \(\ell_1\), \(\ell_2\), or \(\ell_\infty\) norm, and to investigate and extend the limits of the IP approach.

While our work focuses on verifying a given BNN after it has been trained, for completeness, we briefly mention some methods for training BNNs. [1] trained BNNs by using a gradient descent method similar to the training of feedforward neural networks, and [17] trained BNNs by using a pseudo-gradient method. [18] and [19] trained BNNs by solving constraint programming or mixed integer programming problems. To improve training BNNs, [20] applied activation scaling, and [21] applied activation approximation.

Finally, we mention that the BNN verification problem is closely related to the problem of finding counterfactual explanations for such networks. Counterfactual explanations seek to find a feature vector that is as close as possible to a given feature vector while leading to a specific alternative desired prediction ([22], [23]). Recently, counterfactual explanations for several machine learning prediction models, including k-nearest neighbors, have been explored ([24], [25]). A method for solving the BNN verification problem can be used to solve the counterfactual explanation problem by conducting a binary search to (approximately) find the smallest \(\epsilon\) such that the feature is classified as a target class.

This paper is organized as follows. In Section 2, we introduce an existing IP formulation of the BNN verification problem and describe methods for obtaining linear constraints and a linear objective in this IP formulation, including our method for obtaining a linear objective by incorporating the decision of which alternative class a perturbed feature vector is misclassified into a single optimization problem. In Section 3, we describe our technique for deriving layerwise derived valid inequalities and our approach for using these valid inequalities to solve the BNN verification problem. We present the results of a computational study in Section 4 and conclude with a discussion of potential future directions in Section 5.

Notation. We use \([M]\) to represent the set \(\{1, \ldots, M\}\) for \(M \in \mathbb{N}\).

2 IP Formulations↩︎

To obtain an IP formulation for problem 1 , \(\mathbf{x}^0 \in \frac{1}{q} \mathbb{Z}_+^{n^0} \cap [0, 1]^{n^0}\) is defined as a vector of the decision variables that represent a perturbed feature vector, and \(\mathbf{x}^\ell\in \{0, 1\}^{n^\ell}\) is defined as the vector of decision variables that represent the output of the \(\ell^{\mathrm{th}}\) hidden layer for \(\ell\in [L]\), which is obtained recursively from \(\mathbf{x}^{\ell- 1}\) as explained in the following paragraphs. Also, \(X^0\) is defined as the set of all perturbed feature vectors to verify over: \[X^0 = \Big\{\mathbf{x}^0 \in \frac{1}{q} \mathbb{Z}_+^{n^0} \cap [0, 1]^{n^0}: \lVert\mathbf{x}^0 - \bar{\mathbf{x}}\rVert_p \le \epsilon\Big\}.\]

The usual description of a BNN describes the output of each neuron as \(+1\) or \(-1\). In order to derive an IP formulation that can be directly solved by IP solvers, we will instead encode the outputs with binary (\(0/1\)) valued variables. Given a binary output vector \(\mathbf{x}^{\ell-1}\), it can be converted to a \(+1/-1\) valued vector via the transformation \(2 \mathbf{x}^{\ell-1} - \mathbf{1}\). Thus, for each layer \(\ell\in [L+ 1]\), the first component of BNN propagation is to multiply the input \(2 \mathbf{x}^{\ell- 1} - \mathbf{1}\) with \(\mathbf{W}^\ell\) and add \(\mathbf{b}^\ell\). We introduce the notation \(a^\ell= (a_1^\ell, \ldots, a_{n^\ell}^\ell)\) for this affine function of \(\mathbf{x}^{\ell- 1}\): \[a^\ell(\mathbf{x}^{\ell- 1}) = \mathbf{W}^\ell(2 \mathbf{x}^{\ell- 1} - \mathbf{1}) + \mathbf{b}^\ell.\] Then, the output of each each \(\ell\in [L]\) is obtained from its input \(\mathbf{x}^{\ell- 1}\) via the transformation \[\mathbf{x}^\ell= \mathtt{1}_{\mathbb{R}_+}(a^\ell(\mathbf{x}^{\ell- 1})),\] where \(\mathtt{1}_{\mathbb{R}_+}\) is the indicator function mapping non-negative numbers to \(1\) and negative numbers to \(0\).¹

Applying \(a^\ell\) and \(\mathtt{1}_{\mathbb{R}_+}\) alternately \(L\) times and \(a^{L+ 1}\) lastly results in \(f\). Thus, \(f\) can be written as follows: \[f(\mathbf{x}^0) = a^{L+ 1}(\mathtt{1}_{\mathbb{R}_+}(a^L(\mathtt{1}_{\mathbb{R}_+}(a^{L- 1}( \cdots \mathtt{1}_{\mathbb{R}_+}(a^1(\mathbf{x}^0)) \cdots ))))).\]

Using these observations, 1 can be written as the following (nonlinear) IP problem:

2.1 Constraint Formulation↩︎

The layer propagation constraints [eq:origFormPropConstr] and the input perturbation constraint [eq:origFormXInpVar] include nonlinear constraints in their description. We review how these constraints can be reformulated using linear constraints in order to obtain a linear IP formulation.

2.1.1 Layer Propagation Constraint Formulation↩︎

To formulate the layer propagation constraints [eq:origFormPropConstr] with linear constraints, the following constants are defined for \(\ell\in [L]\) and \(i\in N^\ell\): \[\begin{align} &\mathrm{LB}_i^\ell:= \sum_{j\in N^{\ell- 1}} (W_{ij}^\ell- |W_{ij}^\ell|)\\ &\mathrm{UB}_i^\ell:= \sum_{j\in N^{\ell- 1}} (W_{ij}^\ell+ |W_{ij}^\ell|)\\ &R_i^\ell:= \begin{cases} \frac{2}{q} \bigg\lceil\frac{q(\sum_{j\in N^0} W_{ij}^1 - b_i^1)}{2}\bigg\rceil - \frac{1}{q} & (\ell= 1)\\ 2 \bigg\lceil\frac{\sum_{j\in N^{\ell- 1}} W_{ij}^\ell- b_i^\ell}{2}\bigg\rceil - 1 & (\ell\in \{2, \ldots, L\}) .\end{cases}. \end{align} \label{eq:const}\tag{2}\] Since \(W_{ij}^\ell\in \{-1, 0, 1\}\) it follows that for any \(\mathbf{x}^{\ell- 1} \in [0, 1]^{n^{\ell- 1}}\) \[\mathrm{LB}_i^\ell \leq 2 \sum_{j\in N^{\ell- 1}} W_{ij}^\ell x_j^{\ell- 1} \leq \mathrm{UB}_i^\ell.\] The following lemma describes how the layer propagation constraints can be modeled with linear constraints. This extends the formulation of [7], which assumes \(W_{ij}^\ell\in \{-1, 1\}\) and \(b_i^\ell= 0\). The proof is in the Appendix.

Thus, the nonlinear constraints [eq:origFormPropConstr] can be replaced with the linear constraints 3 - 6 .

2.1.2 Input Perturbation Constraint Formulation↩︎

The input perturbation constraint [eq:origFormXInpVar] can be formulated by defining the integer-valued decision variable \(y_i\) to represent \(qx_i^0\), for \(i\in N^0\). Then, \(\mathbf{x}^0 \in X^0\) if and only if \(\mathbf{x}^0 \in \mathbb{R}_+^{n^0}\), and there exist \(\mathbf{y}\in \mathbb{Z}_+^{n^0}\) satisfying the following constraints: \[\begin{align} &\mathbf{y}= q\mathbf{x}^0,\tag{7}\\ &\mathbf{y}\le q\mathbf{1},\tag{8}\\ &\lVert\mathbf{x}^0 - \bar{\mathbf{x}}\rVert_p \le \epsilon.\tag{9} \end{align}\] In the case of \(p= 1\) or \(p= \infty\), 9 can be formulated with linear constraints using standard techniques. In the case of \(p= 2\), 9 can be formulated as a convex quadratic constraint.

2.2 Objective Formulation↩︎

Finally, to obtain a linear IP formulation of [eq:origForm] we explore two methods to obtain a linear objective.

The first approach, which has been used in [7], [8], is to consider separately each alternative class \(t\in N^{L+ 1} \setminus \{\bar{t}\}\). Each such class \(t\) results in the following IP problem with linear objective:

The optimal solution of [eq:origForm] can be obtained by solving [eq:indivForm] for each \(t\in N^{L+ 1} \setminus \{\bar{t}\}\) and choosing \(t\) that achieves the maximum of \(z^*_\epsilon(\bar{\mathbf{x}}, t)\).

For problems with more than two classes, we propose an alternative approach for obtaining a linear objective which directly includes the choice of the alternative class \(t\) in the formulation, thereby avoiding the need to solve [eq:indivForm] for each \(t\in N^{L+ 1} \setminus \{\bar{t}\}\). Thus, for each \(t\in N^{L+ 1} \setminus \{\bar{t}\}\) let \(z_t\in \{0, 1\}\) be a decision variable indicating whether \(t\) is selected as the alternative class, and for each \(i\in N^L\) let \(v_{ti} \in \{0, 1\}\) be a decision variable that is used to represent the product \(z_tx_i^L\). The proof of the following lemma is in the appendix.

2.3 Valid Inequalities for Neurons↩︎

[7] explored the following convex hull for a single neuron for \(\ell\in \{2, \ldots, L\}\) and \(i\in N^\ell\) in the case of \(\mathbf{b}^\ell= \mathbf{0}\): \[\mathrm{conv}(\{(\mathbf{x}^{\ell- 1}, x_i^\ell) \in \{0, 1\}^{n^{\ell- 1}} \times \{0, 1\}: x_i^\ell= \mathtt{1}_{\mathbb{R}_+}(a_i^\ell(\mathbf{x}^{\ell- 1}))\}). \label{eq:singleConvHull}\tag{10}\]

In the following theorem, we extend the results in [7] to a general case where \(\mathbf{b}^\ell\) are not necessarily zero. This theorem is proved in the Appendix.

Since the number of inequalities described in this theorem grows exponentially with the number of nodes in a hidden layer, these would be used to improve the LP relaxation of the formulation [eq:incorpForm] by adding them via a cutting plane procedure.

[8] extended [7] by investigating the following convex hull for a pair of neurons in the same hidden layer for \(\ell\in \{2, \ldots, L\}\) and \(i, k\in N^\ell\) satisfying \(i< k\) in the case where \(n^{\ell- 1}\) is an even integer larger than \(2\), every entry in \(\mathbf{W}^\ell\) is nonzero, and \(\mathbf{b}^\ell= 0\): \[\mathrm{conv}(\{(\mathbf{x}^{\ell- 1}, x_i^\ell, x_k^\ell): \mathbf{x}^{\ell- 1} \in \{0, 1\}^{n^{\ell- 1}}, x_i^\ell= \mathtt{1}_{\mathbb{R}_+}(a_i^\ell(\mathbf{x}^{\ell- 1})), x_k^\ell= \mathtt{1}_{\mathbb{R}_+}(a_k^\ell(\mathbf{x}^{\ell- 1}))\}).\] In this work, a relaxation IP problem to the IP problem [eq:indivForm] is solved for a fixed alternative class \(t\in N^{L+ 1}\) to tackle the BNN verification problem. This computational study shows some improvement in the optimal objective value of the LP relaxation of the IP formulation compared to [7], but this does translate to significant improvement in the ability to solve the BNN verification problem.

3 Layerwise Derived Valid Inequalities↩︎

We explore an alternative approach for deriving valid inequalities for [eq:indivForm] that is based on recursively approximating the set of “reachable vectors” at each layer.

Starting from \(X^0\), \(X^\ell\) is defined recursively as follows for \(\ell\in [L]\): \[X^\ell= \mathtt{1}_{\mathbb{R}_+}(a^\ell(X^{\ell- 1})).\] Then, \(X^\ell\) is the set of all output vectors of layer \(\ell\) that can be obtained from a perturbed input vector in the set \(X^0\). Since the objective [eq:origFormObj] only depends on \(\mathbf{x}^L\), [eq:origForm] can be formulated as the following IP problem:

The above IP problem implies that [eq:origForm] can be solved with access to a description for \(X^L\), which motivates studying valid inequalities for this set. We observe that valid inequalities for \(X^\ell\) can be obtained using a relaxation of the set \(X^{\ell-1}\) since \(X^\ell= \mathtt{1}_{\mathbb{R}_+}(a^\ell(X^{\ell- 1}))\). Since a full description for \(X^0\) is given, this motivates an iterative approach in which we derive valid inequalities for the set \(X^{\ell}\) based on an outer approximation of the set \(X^{\ell-1}\) for \(\ell=2,\ldots,L\). We refer to such inequalities as layerwise derived valid inequalities.

3.1 Variable Fixing↩︎

We first investigate layerwise derived valid inequalities for \(X^\ell\) for \(\ell\in [L]\) of the form \[x_i^\ell= \frac{1-c_i}{2} \label{eq:varFix}\tag{11}\] for \(i\in N^\ell\) and \(c_i\in \{-1, 1\}\). In the case of \(c_i= -1\), 11 is equivalent to \(x_i^\ell= 1\) and 11 is equivalent to \(x_i^\ell= 0\) in the case of \(c_i= 1\). For these reasons, valid equalities for \(X^\ell\) of the form 11 are called variable fixings. This idea of identifying variable fixings is related to the work of [9] who apply a similar technique to fix binary variables in the IP formulation of the feedforward neural network verification problem.

Our approach is to consider each equality of the form 11 and determine if it is a valid equality. The following theorem describes an IP problem that can be solved to determine validity of such an equality.

The details are described in Section 3.3, but to preview this, the idea is to use Theorem 2 to obtain outer approximations of the sets \(X^\ell\) by setting \(X_{\mathrm{out}}^0\) as \(X^0\), and then recursively using the theorem to derive valid inequalities to define valid inequalities to define \(X_{\mathrm{out}}^{\ell}\) based on the previously determined \(X_{\mathrm{out}}^{\ell-1}\) for \(\ell=1,\ldots,L\).

To prove this theorem, we use the following lemma.

Although checking whether a variable fixing is valid requires solving the IP ?? , this IP is much easier to solve than the original verification IP as it considers only a single layer at a time and the constraints defining the outer approximation of \(X^{\ell- 1}\) do not involve any “big-\(M\)” constraint ?? . Furthermore, it can be solved very efficiently for certain structures of \(X_{\mathrm{out}}^{\ell-1}\) as discussed in the following paragraphs.

First, consider the case \(\ell\in \{2, \ldots, L\}\) and assume \(X_{\mathrm{out}}^{\ell- 1}\) is defined only by variable fixings, e.g., as is the case if it is defined recursively only the result of Theorem 2. Then, ?? can be solved by setting \(x_j^{\ell- 1}\) to \(1\) if \(c_iW_{ij}^\ell= 1\) and the variable \(x_{j}^{\ell-1}\) is not fixed in \(X_{\mathrm{out}}^{\ell- 1}\), and setting all other unfixed variables to \(0\).

For \(\ell= 1\), the special form of the set \(X^0\) also allows ?? to be solved efficiently. We describe here how this is done for the case \(p= 1\). The methods for solving ?? for the cases \(p=2\) and \(p=\infty\) are described in Appendix 6.4. For \(p=1\), the following inequality holds because \(\mathbf{x}^0 \in [0, 1]^{n^0}\): \[\label{eq:ubcx} c_i\sum_{j\in N^0} W_{ij}^1 x_j^0 \le \sum_{j\in N^0} \max(c_iW_{ij}^1, 0).\tag{12}\] Inequality 12 is satisfied at equality by \(\hat{\mathbf{x}}\) defined by \(\hat{x}_j=\bar{x}_j\) for \(j\) with \(W_{ij}^1 = 0\), \(\hat{x}_j=1\) for \(j\) with \(W_{ij}^1 = c_i\), and \(\hat{x}_j=0\) for \(j\) with \(W_{ij}^1 = \;-c_i\). If \(\|\hat{\mathbf{x}}- \bar{\mathbf{x}}\|_1 \le \epsilon\), \(\hat{\mathbf{x}}\in X^0\), so the optimal objective value is \(\sum_{j\in N^0} \max(c_iW_{ij}^1, 0)\). Otherwise, by iteratively decreasing \(\hat{x}_j\) by \(\frac{1}{q}\) for \(j\) with \(W_{ij}^1=c_i\) and \(\hat{x}_j > \bar{x}_j\) or increasing \(\hat{x}_j\) by \(\frac{1}{q}\) for \(j\) with \(W_{ij}^1=-c_i\) and \(\hat{x}_j < \bar{x}_j\) we can obtain a solution \(\mathbf{x}^* \in X^0\) satisfying \(\|\mathbf{x}^* - \bar{\mathbf{x}}\|_1 = \frac{1}{q}\lfloor q\epsilon\rfloor\), \(x^*_j=\bar{x}_{j}\) for \(j\) with \(W_{ij}^1 = 0\), \(x^*_j>\bar{x}_{j}\) for \(j\) with \(W_{ij}^1 = c_i\), and \(x^*_j< \bar{x}_{j}\) for \(j\) with \(W_{ij}^1 = -c_i\). The following inequalities are satisfied by every \(\mathbf{x}\in X^0\) and hold at equality for \(\mathbf{x}^*\): \[\begin{align} c_i\sum_{j\in N^0} W_{ij}^1 x_j^0 &= \sum_{j\in N^0} c_iW_{ij}^1 (x_j^0 - \bar{x}_j) + c_i\sum_{j\in N^0} W_{ij}^1 \bar{x}_j&\\ &\le \sum_{j\in N^0} |x_j^0 - \bar{x}_j| + c_i\sum_{j\in N^0} W_{ij}^1 \bar{x}_j&(c_iW_{ij}^1 \in \{-1, 0, 1\} )\\ &\le \frac{1}{q}\lfloor q\epsilon\rfloor + c_i\sum_{j\in N^0} W_{ij}^1 \bar{x}_j& \end{align}\] where the last inequality follows because \(\sum_{j\in N^0} |x_j^0 - \bar{x}_j| = \|\mathbf{x}^0 - \bar{\mathbf{x}}\|_1 \le \epsilon\) and \(\sum_{j\in N^0} |x_j^0 - \bar{x}_j|\) is a multiple of \(\frac{1}{q}\). It follows that the optimal objective value of ?? in this case is \(\frac{1}{q}\lfloor q\epsilon\rfloor + c_i\sum_{j\in N^0} W_{ij}^1 \bar{x}_j\) and this value is achieved by \(\mathbf{x}^*\).

3.2 Two-variable Inequalities↩︎

We next seek to derive sufficient conditions for inequalities of the following form to be valid inequalities for \(X^\ell\) for \(\ell\in [L]\): \[c_ix_i^\ell+ c_kx_k^\ell\le \frac{c_i+ c_k}{2}, \label{eq:twoVarIneq}\tag{13}\] for \(i, k\in N^\ell\) satisfying \(i> k\), and \(c_i, c_k\in \{-1, 1\}\). The following theorem presents an IP that yields a sufficient condition for 13 to be a valid inequality for \(X^\ell\).

Once again, solving 13 is expected to be much simpler than solving the original verification problem as it only includes the constraints defining the outer approximation of \(X^{\ell- 1}\) and one other constraint and does not contain any “big-\(M\)” constraints.

3.3 Algorithms↩︎

We now describe how we propose to use the layerwise derived valid inequalities to solve the IP problem for BNN verification.

Our proposed method is described in Algorithm 1. The algorithm consists of two phases. In the first phase, we only try to fix variables, since in this case the validity verification problem ?? is easy to solve. This is described in lines [p1start] - [p1end], where variable fixings are identified by calling Algorithm 2 for each \(\ell\in [L]\), which outputs an outer approximation for \(X^\ell\) (\(X_{\mathrm{out}}^\ell\)) based on fixing variables in the set \(N_{\mathrm{fix}}^\ell\subseteq N^\ell\). Then, the generated variable fixings are added to the IP formulation [eq:incorpForm] as constraints and [eq:incorpForm] is solved by an IP solver with a limit of one node (i.e., the IP solver only processes the root node of the branch-and-bound tree). If the IP solver solves the verification problem within that limit (e.g., by finding a feasible solution with positive objective value or proving a non-positive upper bound), then it returns the result accordingly. Otherwise, Algorithm 1 moves to the second phase. In this phase, two-variable inequalities are generated by Algorithm 4 for \(\ell\in [L]\) and additional variable fixings are checked, which may be possible due to the improved outer approximations defined by the two-variable inequalities. For layers \(\ell\in \{2,\ldots,L\}\), Algorithm 4 outputs an outer approximation \(X_{\mathrm{out}}^\ell\) that is potentially a subset of the one generated in the first phase. The inequalities defining these outer approximations are added to the IP formulation [eq:incorpForm] and it is solved. In order to answer the verification question, it is only necessary to determine the sign of the optimal value of [eq:incorpForm]. Thus, when solving [eq:incorpForm] the process can be terminated if either the best objective value of a feasible solution found is positive or the best upper bound for the optimal objective value is non-positive.

Algorithms 2, 3, and 4 generate variable fixings in phase 1, variable fixings in phase 2, and two-variable inequalities, respectively. These algorithms check whether each candidate layerwise derived valid inequality is valid for \(X^{\ell}\) by solving the IP subproblem ?? or 13 . When solving these IP subproblems, if the upper bound on the optimal objective value falls below \(c_iR_i^\ell/2\) the optimization process is terminated since we then know the candidate inequality is valid. On the other hand, if the lower bound (incumbent objective value) exceeds \(c_iR_i^\ell/2\), the optimization process is terminated since we can then conclude that we are not able to verify validity of the candidate inequality. In Algorithms 2 and 3, the neuron fixed by the candidate is added to \(N_{\mathrm{fix}}^\ell\) because this set is exploited in Algorithms 3 and 4 to retrieve candidates for layerwise derived valid inequalities.

To limit the time spent checking validity of candidates, Algorithms 3 and 4 use an inner approximation \(X_{\mathrm{in}}^\ell\) for \(\mathtt{1}_{\mathbb{R}_+}(a^\ell(X_{\mathrm{out}}^{\ell- 1}))\) to detect candidates that will not lead to a valid inequality, and hence avoid solving ?? or 13 for those candidates. For each candidate layerwise derived inequality, if the maximum of the left-hand side over \(X_{\mathrm{in}}^\ell\) is larger than the right-hand side, this candidate is not valid for \(X^{\ell}\), so solving the IP subproblem can be avoided. In Algorithm 1, \(X_{\mathrm{in}}= (X_{\mathrm{in}}^0, \cdots, X_{\mathrm{in}}^L)\) is initialized by propagating the input feature vector \(\bar{\mathbf{x}}\) through the BNN. In Algorithm 2 these sets are expanded by propagating the optimal solution obtained when solving ?? at layer 1 (which are in \(X^0\) by definition) through the BNN, thus yielding solutions in \(X^\ell\) which are added to \(X_{\mathrm{in}}^\ell\) for \(\ell=0,\ldots,L\). In Algorithms 3 and 4, for each layer \(\ell\in \{1,\ldots, L\}\), each time we solve ?? or 13 to check validity of a candidate inequality we collect the set of feasible solutions \(X_{\mathrm{feas}}^{\ell-1}\) obtained by the solver, which by definition are in the set \(X_{\mathrm{out}}^{\ell-1}\). These are then propagated forwarded as \(X_{\mathrm{feas}}^{\ell'} = \mathtt{1}_{\mathbb{R}_+}(a^{\ell'}(X_{\mathrm{feas}}^{\ell' - 1}))\) for \(\ell'=\ell,\ldots,L\) and these sets are added to \(X_{\mathrm{in}}^{\ell'}\). Since the set \(X_{\mathrm{out}}^{\ell-1}\) is completely determined when processing layer \(\ell\) in Algorithm 3, this process ensures that \(X_{\mathrm{in}}^{\ell} \subseteq X_{\mathrm{out}}^{\ell}\) for all \(\ell\), and hence if a candidate valid inequality is violated by a vector in \(X_{\mathrm{in}}^{\ell}\) we can conclude that solving ?? or 13 cannot lead to a verification that the inequality is valid.

As a final strategy to limit the computational time spent solving 13 on candidates that do not yield valid inequalities, Algorithm 4 includes a heuristic stopping rule that quits the process if it seems unlikely to yield more additional valid inequalities. The idea is to check the candidate inequalities in a sequence defined by a score that correlates with how likely they are to yield a valid inequality, and then terminate once the number of consecutive failed attempts exceeds a pre-specified limit ‘\(\mathrm{max\_fail}\)’. The score that we use for a candidate inequality \((i, c_i, k, c_k) \in N_{\mathrm{pair}}^\ell\) also leverages the inner approximations we have built and is defined as: \[\text{score}(i, c_i, k, c_k) = \frac{|\{\hat{\mathbf{x}}^\ell\in X_{\mathrm{in}}^\ell: \hat{x}_i^\ell= \frac{-c_i+ 1}{2}\}|+ |\{\hat{\mathbf{x}}^\ell\in X_{\mathrm{in}}^\ell: \hat{x}_k^\ell= \frac{-c_k+ 1}{2}\}|}{|X_{\mathrm{in}}^\ell|}.\] To understand the intuition for this score, consider the case of \(c_i=c_k=1\), so that the candidate inequality is of the form \(x^\ell_i + x^\ell_k \leq 1\), and the score sums the fraction of solutions in the inner approximation \(X_{\mathrm{in}}^\ell\) which have \(x^\ell_i = 0\) and which have \(x^\ell_k = 0\). A high score suggests most solutions in \(X_{\mathrm{out}}^\ell\) have either \(x^\ell_i = 0\) or \(x^\ell_k = 0\) and hence the candidate inequality might be satisfied. The intuition behind the other cases is similar.

4 Computational Study↩︎

We pursue a computational study to investigate how IP methods work to solve the BNN verification problem for multiple test instances. The following five IP methods are considered in the computational study:

• Many-IP: solve the IP problem [eq:indivForm] for all \(t\),

• 1-IP: solve the IP problem [eq:incorpForm],

• 1-IP+HG: solve [eq:incorpForm] by employing a constraint generation approach with ?? ,

• 1-IP+Fix: solve [eq:incorpForm] by employing the variant of Algorithm 1 where non-root nodes are explored in solving [eq:incorpForm] in phase 1, and phase 2 is skipped,

• 1-IP+Fix+2Var: solve [eq:incorpForm] by employing Algorithm 1.

Many-IP and 1-IP are compared to explore the impact of the technique for obtaining a linear objective. The methods 1-IP, 1-IP+HG, 1-IP+Fix, and 1-IP+Fix+2Var are compared to explore the impact of the techniques for generating layerwise derived valid inequalities.

4.1 Test Instances↩︎

Unless stated otherwise, the maximum perturbation from \(\bar{\mathbf{x}}\) is defined using an \(\ell_1\) norm in our test instances.

Each test instance in the computational study consists of a BNN, \(\bar{\mathbf{x}}\), and \(\epsilon\). Nine BNNs are trained for the computational study by using [1]’s method based on a gradient descent method with the MNIST train dataset. The MNIST dataset consists of feature vectors representing handwritten digits from 0 to 9, so \(n^{L+ 1} = 10\). These feature vectors are originally 784-dimensional non-negative integer vectors whose coordinates are at most 255, but they are scaled to 784-dimensional non-negative real vectors whose coordinates are at most 1, so \(n^0 = 784\) and \(q= 255\). The number of hidden layers, the number of neurons in hidden layers, and the error rate of each network are reported in Table 1. The error rate of each network is computed as the portion of feature vectors misclassified by this network in the MNIST test dataset.

For each network, we consider 10 instances, defined by 10 different feature vectors \(\bar{\mathbf{x}}\). For each digit from 0 to 9, one feature vector in the MNIST test dataset whose \(\bar{t}\) is this digit is randomly chosen.

For each network and feature vector \(\bar{\mathbf{x}}\), six values are obtained for \(\epsilon\) by employing Algorithm 6 with \(\mathrm{max\_iter}= 6\). Algorithm 6 is a binary search that attempts to find the largest \(\epsilon\) under which \(\bar{\mathbf{x}}\) can be \(\epsilon\)-verified using a given method and time limit in each iteration. In Algorithm 6, if \(\bar{\mathbf{x}}\) is \(\epsilon\)-verified for the current \(\epsilon\), \(\epsilon\) is increased based on the smallest input perturbation \(\epsilon_{\mathrm{UB}}\) under which \(\bar{\mathbf{x}}\) may not be \(\epsilon_{\mathrm{UB}}\)-verified. Otherwise, \(\epsilon\) is decreased to the average of \(\epsilon\) and the largest input perturbation \(\epsilon_{\mathrm{LB}}\) under which \(\bar{\mathbf{x}}\) is \(\epsilon_{\mathrm{LB}}\)-verified. Different IP methods for the BNN verification may result in different sequences of \(\epsilon\) in Algorithm 6 based on their ability to solve the verification problem at a given \(\epsilon\) within the time limit. This approach is used for determining the \(\epsilon\) values in the test instances in order to find challenging instances, i.e., those in which determining whether \(\bar{\mathbf{x}}\) is \(\epsilon\)-verified is nontrivial. The initial value \(\epsilon_{\mathrm{init}}\) is set to \(1\) in the case that the maximum perturbation from \(\bar{\mathbf{x}}\) is defined by an \(\ell_1\) norm, \(\frac{1}{255}\) in the case that the maximum perturbation is defined by \(\ell_\infty\) norm, and \(\frac{1}{32}\) in the case that the maximum perturbation is defined by \(\ell_2\) norm.

4.2 Implementation Details↩︎

In the IP method 1-IP+HG, ?? are generated by solving the LP relaxation problem of the IP problem [eq:incorpForm] and adding violated inequalities ?? and ?? in an iterative manner. In each iteration, for each \(\ell\in \{2, \ldots, L\}\) and \(i\in N^\ell\), an inequality ?? with the largest violation is added to [eq:incorpForm]. Likewise, ?? with the largest violation is added to [eq:incorpForm] as a constraint. This iteration is repeated until either no violated inequalities are found or the optimal objective value of the LP relaxation problem does not improve by 1% over the last 10 iterations. Our test instances are different than those used in [7] (ours have more layers in the BNN and non-binary inputs) and we obtain qualitatively different results than those presented in [7]. Thus, to validate our implementation, in Appendix 6.6 we present results of additional experiments with the method 1-IP+HG on the instances used in [7] which indicates that our implementation achieves qualitatively similar results to what is reported in [7] on those instances.

In the IP method 1-IP+Fix and 1-IP+Fix+2Var, in Algorithm 2, the IP subproblem ?? is solved using the method described at the end of Section 3.1 (i.e., without an IP solver).

In 1-IP+Fix+2Var, two-variable inequalities for the first hidden layer are not generated, and variable fixings for the second hidden layer are not generated in phase 2 because two-variable inequalities for the first hidden layer are unlikely to be generated for the test instances. The networks are dense, and \(\epsilon\) is much smaller than \(n^0 = 784\) in test instances, so for \(i, k\in N^1\) satisfying \(i< k\) and \(c_i, c_k\in \{-1, 1\}\), it is likely to find \(\hat{\mathbf{x}}^0 \in X^0\) where \(2 c_i\sum_{j\in N^0} W_{ij}^1 \hat{x}_j^0\) and \(2 c_k\sum_{j\in N^0} W_{kj}^1 \hat{x}_j^0\) achieve their maximum over \(X^0\). If both \(x_i^1\) and \(x_k^1\) are not fixed, \(\hat{\mathbf{x}}^0\) violates 13 , so two-variable inequalities for the first hidden layer are rarely generated. As a result of not generating two-variable inequalities for the first hidden layer, variable fixings for the second hidden layer cannot be generated in phase 2. We use \(\mathrm{max\_fail}= 100\) as the number of consecutive failures after which we terminate attempting to generate two-variable inequalities.

The time limit for all IP methods is 3600 seconds. In the IP method Many-IP, the time limit to solve the IP problem [eq:indivForm] for each \(t\) is set to 400 seconds because there are nine alternative classes. In 1-IP+HG, 1-IP+Fix, and 1-IP+Fix+2Var, a time limit on the phase for generating valid inequalities of 2700 seconds is imposed. The time limit to solve [eq:incorpForm] after generating the inequalities is set to 3600 seconds subtract the time spent generating valid inequalities. If the time limit is hit, the best upper bound for the optimal objective value is used instead of the best objective value to answer the BNN verification problem because \(\bar{\mathbf{x}}\) is \(\epsilon\)-verified if and only if the best upper bound is non-positive.

All IP methods are implemented in Python, and Gurobi 10.0.3 is used as the IP solver. All experiments in the computational study are run on an Ubuntu desktop with 32 GB RAM and 16 Intel Core i7-10700 CPUs running at 2.90 GHz.

4.3 Results↩︎

We first collect the maximum \(\epsilon\) under which \(\bar{\mathbf{x}}\) is \(\epsilon\)-verified for each IP method, network, and \(\bar{\mathbf{x}}\) in the case that the maximal perturbation from \(\bar{\mathbf{x}}\) is defined using an \(\ell_1\) norm. The maximum \(\epsilon\) is obtained from using Algorithm 6 with \(\mathrm{max\_iter}= 6\). For the same test instance, one IP method may succeed in verifying the BNN, but another IP method fails because it may not be able to verify the BNN within the time limit. As a result, the set of \(\epsilon\) defining the six instances may differ by IP methods for the same network and \(\bar{\mathbf{x}}\), and the maximum \(\epsilon\) may differ.

For each IP method and network, the arithmetic mean of the maximum \(\epsilon\) over ten \(\bar{\mathbf{x}}\) is presented in Table 2. Many-IP achieves the smallest maximum \(\epsilon\), and 1-IP+HG achieves a smaller maximum \(\epsilon\) than the 1-IP. Except for Network 4, 1-IP+Fix results in a larger maximum \(\epsilon\) than 1-IP. Method 1-IP+Fix+2Var yields the largest maximum \(\epsilon\) in every case. These results indicate that the using the proposed variable fixing approach yields modest improvement in the ability to find the maximum \(\epsilon\) at which an input can be verified, and using two-variable layer-wise derived inequalities yields significantly more improvement.

We next investigate in more detail the ability of the different methods to solve verification problems for various values of \(\epsilon\). We first compare Many-IP and 1-IP. To do so, for each network and feature vector we consider the subset of \(\epsilon\) values that are solved by both of these methods, and compute the shifted geometric mean of solution times to solve the IP problem [eq:incorpForm] (for 1-IP) or to solve [eq:indivForm] for all \(t\) (for Many-IP), with a shift of 1. These results are presented in Table 3, where for Many-IP, the number in parenthesis is the number of test instances where the time limit is hit in solving [eq:indivForm] for at least one \(t\) and for 1-IP, the number in parenthesis is the number of test instances where the time limit is hit in solving [eq:incorpForm]. We find that 1-IP solves the verification problem much more quickly than Many-IP, indicating that our technique for obtaining a linear objective leads to faster BNN verification. Considering the superiority of 1-IP over Many-IP, we exclude Many-IP from further comparisons.

In Table 4 we report the geometric average verification times over the instances tested by all IP methods based on 1-IP. The set of test instances in this case is based on the set of \(\epsilon\) values that were solved by all 1-IP methods. This set of test instances may differ from the set used in the Many-IP and 1-IP comparison due to the exclusion of Many-IP from this comparison. The number in parenthesis is the number of test instances where the time limit is hit when solving [eq:incorpForm]. We find that 1-IP+HG’s verification times are larger than 1-IP’s, and 1-IP+Fix’s are shorter than 1-IP’s. Thus, we conclude that using inequalities ?? does not accelerate BNN verification, but variable fixings do. We also find that 1-IP+Fix+2Var’s verification times are similar to 1-IP+Fix’s for Networks 1-4, and 1-IP+Fix+2Var’s verification times are shorter than 1-IP+Fix’s for Networks 5-9. This suggests that the two-variable layerwise derived inequalities are more helpful for BNN verification of networks with more hidden layers.

We next investigate the LP gaps of [eq:incorpForm] over instances tested by all IP methods based on 1-IP. The LP gap is defined as \((z_{\mathrm{LP}}^*- z^*)/\bar{z}_{\mathrm{UB}}\) where \(z_{\mathrm{LP}}^*\) is the optimal objective value of the LP relaxation problem of [eq:incorpForm] after adding any valid inequalities introduced by the method, \(z^*\) is the best objective value from a feasible solution to [eq:incorpForm] obtained by the method, and \[\bar{z}_{\mathrm{UB}}= \max_{\substack{\mathbf{x}^L\in \{0, 1\}^{n^L}\\ t\in N^{L+ 1}}} (a_t^{L+ 1}(\mathbf{x}^L) - a_{\bar{t}}^{L+ 1}(\mathbf{x}^L)).\] We use \(\bar{z}_{\mathrm{UB}}\) as the denominator instead of \(|z^*|\) because \(z^*\) may be near zero for some instances, which would skew interpretation of these relative optimality gaps when averaged over different instances. The quantity \(\bar{z}_{\mathrm{UB}}\) is positive by definition and is a natural upper bound on the optimal value of [eq:incorpForm]. Table 5 reports the average LP gaps for each method. We find that the LP gaps of 1-IP+HG and 1-IP are similar, 1-IP+Fix yields significant improvement over 1-IP, and 1-IP-Fix+2Var yields further significant improvement beyond that. Although 1-IP+HG does not improve LP gaps on these instances, we remark that in the supplemental experiments reported in Appendix 6.6 we found that on the instances used in [7] the use of the single-neuron valid inequalities did yield modest reduction in the LP gap, although even on these instances this reduction still did not translate into reduced time to solve the verification problem.

Table 6 presents additional results for the IP methods based on 1-IP. These results are averaged over all 445 instances solved by all these methods. The row ‘Preprocessing time’ presents the shifted geometric mean of times to generate layerwise derived valid inequalities (resp. ?? ) for 1-IP+Fix and 1-IP+Fix+2Var (resp. 1-IP+HG). The row ‘verification time’ is the shifted geometric mean of times to solve the verification problem – if a method does not solve an instance within the time limit, then the time limit is used for that instance. The row ‘# of solved instances (veri.)’ reports the number of instances for which each IP method solved the verification problem within the time limit. Row ‘# of nodes (veri.)’ is the shifted geometric mean of the number of nodes in the branch-and-bound tree when solving the verification problem (again, for instances hitting the time limit, this is just the number of nodes explored within the limit). In all the verification times the solution of [eq:incorpForm] is terminated as soon as the verification question is answered (e.g., if the upper bound becomes non-positive or the lower bound becomes positive). To provide further insight into the quality of these formulations we also attempted to solve [eq:incorpForm] without terminating once the verification question is answered. The results of this experiment are reported in the rows ‘Optimization time’, ’ (# of solved instances (opti.)‘,’(# of nodes (opti.))‘, which are defined analogously as the verification results. Row ’optimality gap’ presents the arithmetic mean of relative optimality gaps obtained in solving [eq:incorpForm] without the early termination, where the relative optimality gap is defined as \((\bar{z}- z^*)/\bar{z}_{\mathrm{UB}}\), where \(\bar{z}\) is the best bound for the optimal objective value of [eq:incorpForm] obtained by the method. A shift is set to 1 in every shifted geometric mean. These results confirm the significant time reduction obtained using variable fixing and the proposed two-variable inequalities, and indicate that this reduction is obtained thanks to the dramatic reduction in the number of branch-and-bound nodes required to either answer the verification problem or to solve the optimization problem. The results also reinforce that the single-neuron inequalities used in 1-IP+HG do not yield improvement.

All results so far have been reported for the case in which the maximum perturbation from \(\bar{\mathbf{x}}\) is defined using the \(\ell_1\) norm. To illustrate the versatility of the IP approach, we also conducted analogous experiments using the \(\ell_\infty\) and \(\ell_2\) norm. The detailed results of these experiments are presented in Appendix 6.5, and we present a summary here. Table 7 presents the arithmetic mean of the maximum \(\epsilon\) over ten feature vectors \(\bar{\mathbf{x}}\) for each IP method on Network 5. Once gain, Many-IP achieves the smallest maximum \(\epsilon\) in both cases. Method 1-IP+HG again yields smaller maximum \(\epsilon\) than 1-IP in the case of \(p= \infty\) and the same in the case of \(p= 2\). 1-IP+Fix yields similar maximum \(\epsilon\) to 1-IP in the case of \(p= \infty\) and larger maximum \(\epsilon\) than 1-IP in the case of \(p= 2\). Thus, it seems that generating variable fixings alone does not improve BNN verification by the IP methods in the case of \(p= \infty\), but it does in the case of \(p= 2\). 1-IP+Fix+2Var yields the largest maximum \(\epsilon\) in both cases, showing that generating layerwise derived valid inequalities enables verifying the BNN against a higher range of \(\epsilon\) in the case of \(p= \infty\) and \(p= 2\).

5 Conclusion↩︎

In this paper, we investigate an IP method that exploits the technique for obtaining a linear objective and a technique for generating layerwise derived valid inequalities to solve the BNN verification problem. Our computational study shows that our IP method verifies BNNs against a higher range of input perturbation than existing IP methods. The technique for obtaining a linear objective leads to solving the BNN verification problem faster by solving a single IP problem instead of multiple IP problems for each alternative class. The technique for generating layerwise derived valid inequalities enables more efficient BNN verification by yielding smaller LP gaps at the expense of solving IP subproblems involving a single layer.

One future direction on IP methods for the BNN verification problem is investigating valid inequalities that involve more than two variables, as the inequalities we investigated contain at most two variables. Another future direction related to the BNN verification problem is counterfactual explanations for BNNs. As discussed in Section 1, these counterfactual explanations can be obtained by solving the BNN verification problem for various input perturbations. Thus, a possible direction of research for counterfactual explanations for BNNs is to consider whether information from these related instances can be shared to accelerate the overall process. For example, it may be possible to reduce the number of candidates for layerwise derived valid inequalities to check.

6 Appendix↩︎

6.1 Proof of Lemma [lem:linPropConstr]↩︎

We prove following Lemma [lem:linPropConstr] on the formulation for the layer propagation constraint [eq:origFormPropConstr] introduced in Section 2.1.1 with the definition on \(\mathrm{LB}_i^\ell\), \(\mathrm{UB}_i^\ell\), and \(R_i^\ell\) for \(\ell\in [L]\) and \(i\in N^\ell\).

6.2 Proof of Lemma 1↩︎

6.3 Proof of Convex Hull Characterization for a Single Neuron↩︎

We prove Theorem 1 on the characterization for the convex hull 10 for a single neuron introduced in Section 2.3, which is an extension of the main result in [7]. To simplify the statement of this theorem, \(\ell\in \{2, \ldots, L\}\) and \(i\in N^\ell\) are fixed, and the following notations are defined for \(j\in N^{\ell- 1}\) with \(\mathrm{LB}_i^\ell\), \(\mathrm{UB}_i^\ell\), and \(R_i^\ell\) defined in Section 2.1.1: \[\begin{align} \mathbf{x}& := \mathbf{x}^{\ell- 1},\\ u& := x_i^\ell,\\ w_j& := W_{ij}^\ell,\\ \mathrm{LB}& := \mathrm{LB}_i^\ell= \sum_{j\in N^{\ell- 1}} (W_{ij}^\ell- |W_{ij}^\ell|),\\ \mathrm{UB}&:= \mathrm{UB}_i^\ell= \sum_{j\in N^{\ell- 1}} (W_{ij}^\ell+ |W_{ij}^\ell|),\\ R&:= R_i^\ell= 2 \bigg\lceil\frac{\sum_{j\in N^{\ell- 1}} W_{ij}^\ell- b_i^\ell}{2}\bigg\rceil - 1. \end{align}\] Then, this theorem can be restated as follows:

Theorem 1 (Restated) The set of \((\mathbf{x}, u) \in [0, 1]^{n^{\ell- 1}} \times [0, 1]\) satisfying the following inequalities for all \(J\subset \{j\in N^{\ell- 1}: w_j\ne 0\}\) is 10 : \[\label{eq:simpSingleConvHullIneq} \begin{align} &\sum_{j\in J} (w_j(2 x_j- 1) - (2 u- 1)) \ge (R- \mathrm{UB}+ 1)u,\label{eq:simpSingleConvHullIneqLowerBound}\\ &\sum_{j\in J} (w_j(2 x_j- 1) - (2 u- 1)) \le (R- \mathrm{LB}- 1)(-u+ 1).\label{eq:simpSingleConvHullIneqUpperBound} \end{align}\] {#eq: sublabel=eq:eq:simpSingleConvHullIneq,eq:eq:simpSingleConvHullIneqLowerBound,eq:eq:simpSingleConvHullIneqUpperBound}

6.4 Implementation Details for Infinity-norm and 2-norm↩︎

We explain implementation details for the case that the maximum perturbation from \(\bar{\mathbf{x}}\) is defined using an \(\ell_\infty\) norm or \(\ell_2\) norm. These details are about how to solve the IP subproblem ?? without an IP solver in the case of \(\ell= 1\) and how to implement the input perturbation constraint 9 .

We first consider solving ?? for the case of layer \(\ell= 1\) and \(p= \infty\). In this case, the optimal solution is obtained by setting \(x_j^0\) to \(\min(\bar{x}_j+ \epsilon, 1)\) for \(j\) satisfying \(c_iW_{ij}^1 = 1\), setting \(x_j^0\) to \(\max(\bar{x}_j- \epsilon, 0)\) for \(j\) satisfying \(c_iW_{ij}^1 = -1\), and setting \(x_j^0\) to \(\bar{x}_j\) for \(j\) satisfying \(c_iW_{ij}^1 = 0\). All constraints defining \(X^0\) are of the form \(|x_j^0 - \bar{x}_j| \le \epsilon\), so this solution is feasible, and it immediate that no other solution can yield a better objective value.

We next consider solving ?? for the case of layer \(\ell= 1\) and \(p= 2\). First, observe that given a feasible solution \(\hat{\mathbf{x}}^0 \in \frac{1}{q} \mathbb{Z}_+^{n^0} \cap [0, 1]^{n^0}\) to ?? , if for some \(j \in N^0\) it holds that \(c_iW_{ij}^1(\hat{x}_j^0 - \bar{x}_j) \le 0\), then another feasible solution with the same or better objective value can be obtained by replacing \(\hat{x}_j^0\) with \(\bar{x}_j\). Hence, we can restrict our search for an optimal solution of ?? to solutions \(\hat{\mathbf{x}}^0\) that satisfy: \[\begin{align} \hat{x}_j^0 \ge \bar{x}_j, &\quad \text{for j such that } c_iW_{ij}^1 = 1,\\ \hat{x}_j^0 \le \bar{x}_j, &\quad \text{for j such that } c_iW_{ij}^1 = -1,\\ \hat{x}_j^0 = \bar{x}_j, &\quad \text{for j such that } W_{ij}^1 = 0. \end{align}\] Let \(N_0^+ = \{j \in N_0 : W_{ij} \neq 0\}\). Thus, we can restrict our attention to solutions defined by a vector \(z \in \{0,1,\ldots,q\}^{N_0^+}\) as follows: \[\hat{x}_j^{0}(z) = \begin{cases} \min(\bar{x}_j+ \frac{z_j}{q}, 1) & \text{if } c_iW_{ij}^1 = 1, \\ \max(\bar{x}_j- \frac{z_j}{q}, 0) & \text{if } c_iW_{ij}^1 = -1, \\ \bar{x}_j& \text{if } W_{ij}^1 = 0 . \end{cases}\] Observe that the objective in ?? and the expression \(\| \hat{\mathbf{x}}(z) - \bar{\mathbf{x}}\|_2\) are monotone increasing in \(z\), and hence an optimal solution will be such that increasing \(z_j\) for any \(j \in N_0^+\) would be infeasible. Next, for \(z \in \{0,1,\ldots,q\}^{N_0^+}\), define the set \[N_0^S(z) = \Bigl\{ j \in N_0^+ : (c_iW_{ij}^1 = 1 \text{ and } \hat{x}_j^0(z) \leq 1 - 1/q ) \text{ or } (c_iW_{ij}^1 = -1 \text{ and } \hat{x}_j^0(z) \geq 1/q ) \Bigr\}.\] If there exists \(z \in \{0,1,\ldots,q\}^{N_0^+}\) in which \(\hat{\mathbf{x}}(z)\) is feasible to ?? and \(N_0^S(z) = \emptyset\), then \(\hat{\mathbf{x}}(z)\) is optimal to ?? since this solution obtains the best possible objective. Otherwise, we claim that there exists an optimal solution \(\hat{\mathbf{x}}(z)\) that satisfies: \[\max\{ z_j : j \in N_0^S(z)\} - \min \{z_j : j \in N_0^S(z) \} \leq 1.\] Consider a solution that violates this condition. Another feasible solution with the same objective value can be obtained by decreasing \(z_{j_1}\) by one and increasing \(z_{j_2}\) by one, where \(j_1 \in \arg \max\{ z_j : j \in N_0^S(z)\}\) and \(j_2 \in \arg \min \{z_j : j \in N_0^S(z) \}\). Repeating this process will eventually yield a vector \(z\) which satisfies this condition since there only finitely many elements achieving the max and min in the condition, and an element will necessarily be removed from one or the other as long as the max is at least two larger than the min.

These arguments imply that an optimal solution to ?? can be obtained by first finding the maximum \(m\) such that the solution \(\hat{\mathbf{x}}(z)\) is satisfies \(\| \hat{\mathbf{x}}(z) - \bar{\mathbf{x}}\|_2 \leq \epsilon\), where \(z_j = m\) for \(j \in N_0^+\). Having found this solution (e.g., by binary search), one would then iteratively select some \(j \in N_0^S(z)\) and increase \(z_j\), and repeat (without selecting any \(j\) more than once) until no more can be increased without violating \(\| \hat{\mathbf{x}}(z) - \bar{\mathbf{x}}\|_2 \leq \epsilon\).

Finally, we describe how the constraint 9 is formulated when using a MIP solver to solve the verification problem in the case of \(p=2\). We define a decision variable \(u_j\) to represent \(|x_j^0 - \bar{x}_j|\) for each \(j\in N^0\). The following inequalities are added as constraints because 9 is satisfied if and only if there exists \(\mathbf{u}\in \mathbb{R}_+^{n^0}\) satisfying the following inequalities: \[\begin{align} u_j&\ge x_j^0 - \bar{x}_j, \quad \forall j\in N^0,\\ u_j&\ge -x_j^0 + \bar{x}_j, \quad \forall j\in N^0,\\ \sum_{j\in N^0} u_j^2 &\le \epsilon^2. \end{align}\] We remark that we found through preliminary empirical study that this formulation technique was computationally superior to the more direct approach of simply formulating the constraint as follows: \[\sum_{j \in N_0} (x_j^0 - \bar{x}_j)^2 \leq \epsilon^2.\]

6.5 Detailed Computational Results for Infinity-norm and 2-norm↩︎

We present more detailed computational results for the case that the maximum perturbation from \(\bar{\mathbf{x}}\) is defined using an \(\ell_\infty\) norm or \(\ell_2\) norm. Test instances for these cases are formed with Network 5 in Table 1, \(\bar{\mathbf{x}}\), and \(\epsilon\) that are chosen as explained in Section 4.1.

We first investigate the time required to use the IP problem [eq:incorpForm] (resp. [eq:indivForm] for all \(t\)) to solve the the verification problem. The shifted geometric mean (with a shift of 1) of times over all instances tested by both IP method Many-IP and 1-IP is computed for each method and displayed in Table 8 for the case of \(p= \infty\) and \(p= 2\), along with the number of test instances. For Many-IP, the number in parenthesis is the number of test instances where the time limit is hit in solving [eq:indivForm] for at least one \(t\). For 1-IP, the number in parenthesis is the number of test instances where the time limit is hit in solving [eq:incorpForm]. We find that 1-IP’s verification times are smaller than Many-IP’s, indicating that our technique for obtaining a linear objective also leads to faster BNN verification in the case of \(p= \infty\) and \(p= 2\).

Table 9 presents results analogous to those of Table 6 in the main sections, but on instances with \(p= \infty\). We find that in this case, neither 1-IP+HG nor 1-IP+Fix yield improvement compared to 1-IP, so generating ?? and variable fixings do not improve BNN verification. Generating variable fixings decreases LP gaps, but this decrease does not lead to more effective BNN verification in the case of \(p= \infty\). On the other hand, 1-IP+Fix+2Var yields improvement compared to the other IP methods based on 1-IP. It implies that our technique for generating layerwise derived valid inequalities, mainly two-variable inequalities, also leads to more efficient BNN verification in the case of \(p= \infty\).

In the case of \(p= 2\), the same metrics are computed over all 22 instances tested by all IP methods based on 1-IP to compare these IP methods. Table 10 shows these metrics. 1-IP+HG does not result in improvement compared to 1-IP, which means that generating ?? also does not improve BNN verification in the case of \(p= 2\). 1-IP+Fix and 1-IP+Fix+2Var yield similar improvements compared to 1-IP, which implies that generating variable fixings improves BNN verification in the case of \(p= 2\), but on these instances there does not appear to be further improvement from two-variable inequalities.

To better compare 1-IP+Fix and 1-IP+Fix+2Var, we report in Table 11 the results comparing these instances on all 49 instances solved by these two methods (many of these were excluded in Table 10 since they were not solved by 1-IP or 1-IP+HG). From this table we observe that spending more time in generating more layerwise derived valid inequalities, mainly two-variable inequalities, results in reducing LP gaps and solving the verification problem more quickly on these instances.

6.6 Computational Study on Convex Hull Characterization for a Single Neuron↩︎

We pursue an additional computational study to validate our implementation for IP methods exploiting Theorem 1 on convex hull characterization for a single neuron. As we see in in Table 5, the IP method 1-IP+HG employing the valid inequalities presented Theorem 1 does not yield significant decreases in LP gaps compared to the IP method 1-IP. However, using these valid inequalities resulted in smaller LP gaps in [7]. One possible reason for this discrepancy is that ?? cannot be used for \(\ell= 1\) in our computational study because \(q= 255\), but these valid inequalities can be used in the computational study of [7] because \(q= 1\). To examine whether this guess is correct, an additional computational study is conducted on test instances based on ones in [7] to assess the impact of exploiting Theorem 1. In this computational study, the IP method Many-IP and Many-IP+HG, which solves the IP problem [eq:indivForm] by employing a constraint generation approach with ?? as 1-IP+HG for each \(t\), are compared because \(t\) is fixed in [7].

Each test instance consists of a BNN, \(\bar{\mathbf{x}}\), and \(\epsilon\). Six BNNs pretrained from [7] are used in this computational study. These networks are trained on the MNIST training dataset by using [1]’s method, which implies \(n^0 = 784\) and \(n^{L+ 1} = 10\). The numbers of hidden layers and the number of neurons in hidden layers of each network are reported in Table 12. Feature vectors in the MNIST dataset are scaled to binary vectors by converting coordinates smaller than 128 to 0 and the other coordinates to 1, which indicates \(q= 1\). Ten feature vectors are selected for \(\bar{\mathbf{x}}\). For each digit from 0 to 9, one feature vector in the MNIST test dataset whose \(\bar{t}\) is this digit is randomly chosen. Positive integers from 1 to 5 are used as \(\epsilon\).

For each network and \(\epsilon\), the arithmetic mean of optimal objective values of the LP relaxation problem of [eq:indivForm] over all \(\bar{\mathbf{x}}\) and \(t\) (LP value) and the shifted geometric mean of times to solve the verification problem for all \(t\) with a shift of 1 over all \(\bar{\mathbf{x}}\) (verification time) are computed for Many-IP and Many-IP+HG to compare these IP methods. Table 13 reports LP values and verification times. The number in parenthesis for verification time is the number of \(\bar{\mathbf{x}}\) where the time limit is hit in solving [eq:indivForm] for at least one \(t\) for each network and \(\epsilon\). Many-IP+HG yields smaller LP values than Many-IP as in [7], so our additional computational study validates our implementation for IP methods employing the valid inequalities in Theorem 1. In particular, the gap in LP values between Many-IP and Many-IP+HG becomes smaller for larger \(L\) and \(\epsilon\) as in [7]. However, Many-IP+HG results in larger verification times than Many-IP, which implies that even though using the inequalities in Theorem 1 improves LP relaxation values, it does not reduce BNN verification time.

References↩︎