BDD4BNN: A BDD-based Quantitative Analysis Framework for Binarized Neural Networks

Verifying and explaining the behavior of neural networks is becoming increasingly important, especially when they are deployed in safety-critical applications. In this paper, we study verification problems for Binarized Neural Networks (BNNs), the 1-bit quantization of general real-numbered neural networks. Our approach is to encode BNNs into Binary Decision Diagrams (BDDs), which is done by exploiting the internal structure of the BNNs. In particular, we translate the input-output relation of blocks in BNNs to cardinality constraints which are then encoded by BDDs. Based on the encoding, we develop a quantitative verification framework for BNNs where precise and comprehensive analysis of BNNs can be performed. We demonstrate the application of our framework by providing quantitative robustness analysis and interpretability for BNNs. We implement a prototype tool BDD4BNN and carry out extensive experiments which confirm the effectiveness and efficiency of our approach.


Introduction
Deep neural networks (DNNs) have achieved human-level performance in several tasks, and are increasingly being incorporated into various application domains such as autonomous driving [4] and medical diagnostics [49]. Modern DNNs usually contain a great many parameters which are typically stored as 32/64-bit floating-point numbers, and require a massive amount of floating-point operations to compute the output for a single input [56]. As a result, it is often challenging to deploy them on resourceconstrained, embedded devices. To mitigate the issue, quantization, which quantizes 32/64-bit floating-points to low bit-width fixed-points (e.g., 4-bits) with little accuracy loss [21], emerges as a promising technique to reduce resource requirements. In particular, binarized neural networks (BNNs) [26] represent the case of 1-bit quantization using the bipolar binaries ±1. BNNs can drastically reduce memory storage and execution time with bit-wise operations, hence substantially improve the time and energy efficiency. BNNs have been demonstrated to achieve high accuracy for a wide variety of applications [33,48,38].
DNNs have been shown to often lack robustness against adversarial samples. Therefore, various formal techniques have been proposed to analyze DNNs, but most of them focus on real-numbered DNNs only. Verification of quantized DNNs has not been thoroughly explored so far, although recent results have highlighted its importance: it was shown that a quantized DNN does not necessarily preserve the properties satisfied by the real-numbered DNN before quantization [12,20]. Indeed, the fixed-point number semantics effectively yields a discrete state space for the verification of quantized DNNs whereas real-numbered DNNs feature a continuous state space. The discrepancy could invalidate current verification techniques for real-numbered DNNs when they are directly applied to quantized counterparts (e.g., both false negative and false positive could occur). Therefore, specialized techniques are required for rigorously verifying quantized DNNs.
Broadly speaking, the existing techniques for quantized DNNs make use of constraint solving which is based on either SAT/SMT or (reduced, ordered) binary decision diagrams (BDDs). A majority of work resorts to SAT/SMT solving. For the 1-bit quantization (i.e., BNNs), typically BNNs are transformed into Boolean formulas where SAT solving is harnessed [42,10,32,41]. Some recent work also studies variants of BNNs [44,27], for instance, three-valued BNNs. For quantized DNNs with multiple bits (i.e., fixed-points), it is natural to encode them as quantifier-free SMT formulas, e.g., using bit-vector and fixed-point theories [7,20,23], so that off-the-shelf SMT solvers can be leveraged. In another direction, BDD-based approaches currently can tackle BNNs only [50]. In a nutshell, they encode a BNN and an input region as a BDD, based on which various analyses can be performed via queries on the BDD. The crux of the approach is how to generate the BDD efficiently. In the work [50], the BDD is constructed by BDD learning [40], thus, currently limited to toy BNNs (e.g., 64 input size, 5 hidden neurons, and 2 output size) with relatively small input regions.
On the other hand, existing work mostly focuses on qualitative verification, which asks whether there exists an input x (in a specified region) for a neural network such that a property (e.g., local robustness) is violated. In many practical applications, checking only the existence is not sufficient. Indeed, for local robustness, such an (adversarial) input almost surely exists which makes a qualitative answer less meaningful. Instead, quantitative verification, which asks how often a property φ is satisfied or violated, is far more useful as it could provide a probabilistic guarantee of the behavior of neural networks. Such a quantitative guarantee is essential to certify, for instance, certain implementations of neural network based perceptual components against safety standards of autonomous vehicles [28,31]. Quantitative analysis of general neural networks, however, is challenging, hence received little attention and for which the results are rather limited so far. DeepSRGR [68] presented an abstract interpretation based quantitative robustness verification approach for DNNs which is sound but incomplete. For BNNs, approximate SAT model-counting solvers ( SAT) are leveraged [6,43] based on the SAT encoding for the qualitative counterpart. Though probably approximately correct (PAC) style guarantees can be provided, verification cost is usually prohibitively high to achieve higher precision and confidence.
Main contributions. We propose a BDD-based framework BDD4BNN to support quantitative analysis of BNNs. The main challenge is how to efficiently build BDDs from BNNs [43]. In contrast to previous work [50] which is learning-based and largely treats the BNN as a blackbox, we directly encode a BNN and the associated input region into BDDs. In a nutshell, a BNN is a sequential composition of multiple internal blocks and one output block. Each block comprises a handful of layers and captures a function f : {+1, −1} n → {+1, −1} m , where n (resp. m) denotes the number of inputs (resp. outputs) of the block. Technically, the function f can be alternatively rewritten as a function over the standard Boolean domain, i.e., f : {0, 1} n → {0, 1} m . A key stepping-stone of our encoding is the observation that the i-th output y i of the block can be captured by a cardinality constraint of the form n j=1 j ≥ k such that y i ⇔ n j=1 j ≥ k, where each literal j is either x j or ¬x j for the input variable x j , and k is a constant. We then present an algorithm to encode a cardinality constraint n j=1 j ≥ k as a BDD with O((n − k) · k) nodes in O((n − k) · k) time. As a result, the input-output relation of each block can be encoded as a BDD, the composition of which yields the BDD for the entire BNN. A distinguished advantage of our BDD encoding lies in its support of incremental encoding. In particular, when different input regions are of interest, there is no need to construct the BDD of the entire BNN from scratch.
Encoding BNNs as BDDs enables a wide variety of applications in security analysis and decision explanation of BNNs. In this paper, we highlight two of them within our framework, i.e., robustness analysis and interpretability. It was shown that DNNs have been suffering from poor robustness to adversarial examples [55, 46,45]. We consider two quantitative variants of the problem: (1) how many adversarial examples does the BNN have in the input region, and (2) how many of them are misclassified to each class? We further provide an algorithm to incrementally compute the (locally) maximal Hamming distance within which the BNN satisfies the desired robustness properties.
Interpretability is an issue arising as a result of the blackbox nature of DNNs [24,39]. In application domains such as medical diagnosis, understanding the decisions made by DNNs has become a pressing need. We consider two problems: (1) why some inputs are (mis)classified into a class by the BNN and (2) are there any essential features in the input region that are common for all samples classified into a class?
Experimental Results. We implement our framework as a prototype tool BDD4BNN using the CUDD package [54], which scales to BNNs with up to 4 internal blocks, 200 hidden neurons, and 784 input size. To the best of our knowledge, it is the first work to precisely analyze such large BNNs that go significantly beyond the state-of-the-art. The experimental results show that BDD4BNN is significantly more efficient and scalable than the learning-based technique [50]. Furthermore, we demonstrate how BDD4BNN can be used in quantitative robustness analysis and decision explanation of BNNs. For quantitative robustness analysis, our experimental results show that BDD4BNN is considerably (5× to 1, 340×) faster and more accurate than the state-of-the-art approximate SAT-based approach [6]. It can also compute precisely the distribution of predicated classes of the images in the input region as well as the locally maximal Hamming distances on several BNNs. For decision explanation, we show the effectiveness of BDD4BNN in computing prime-implicant explanations and essential features of the given input region for some target classes.
In general, our main contributions can be summarized as follows.
-We introduce a novel algorithmic approach for encoding BNNs into BDDs that exactly preserves the semantics of BNNs, which supports incremental encoding. -We propose a framework for quantitative verification of BNNs and in particular, we demonstrate the robustness analysis and interpretability of BNNs. -We implement the framework as an end-to-end tool BDD4BNN and conduct thorough experiments on various BNNs, demonstrating the efficiency and effectiveness of BDD4BNN. Outline. The remainder of this paper is organized as follows. Section 2 briefly introduces BNNs and BDDs. Section 3 and Section 4 present our BDD-based quantitative analysis framework and its applications respectively. Section 5 reports the evaluation results. Section 6 discusses related work. Finally, we conclude this work in Section 7.

Preliminaries
In this section, we briefly introduce binarized neural networks (BNNs) and (reduced, ordered) binary decision diagrams (BDDs). We denote by R, N, B, and B ±1 the set of real numbers, the set of natural numbers, the standard Boolean domain {0, 1} and the integer set {+1, −1}. For n ∈ N, we denote by [n] the set {1, · · · , n}. We will use W, W , . . . to denote (2-dimensional) matrices, x, v, · · · to denote (row) vectors, and x, v, . . . to denote scalars. We denote by W i,: and W :, j the i-th row and j-th column of the matrix W. Similarly, we denote by x j and W i, j the j-th entry of x and W i,: respectively. In this work, Boolean values 1/0 will be used as integers 1/0 in arithmetic computations without typecasting.

Binarized Neural Networks
A binarized neural network (BNN) [26] is a neural network where weights and activations are predominantly binarized over the domain B ±1 . In this work, we consider feed-forward BNNs. As shown in Figure 1, a BNN can be seen as a sequential composition of several internal blocks and one output block. Each internal block comprises 3 layers: a linear layer (LIN), a batch normalization layer (BN), and a binarization layer (BIN). The output block comprises a linear layer and an ARGMAX layer. Note that the input/output of internal blocks and the input of the output block are all vectors over B ±1 .
±1 → B s with s classes is given by a tuple of blocks ±1 is an internal block comprising a LIN layer t lin i , a BN layer t bn i and a BIN t bin ±1 → B s is the output block comprising a LIN layer t lin d+1 and an ARGMAX layer t am d+1 with t d+1 = t am d+1 • t lin d+1 , Weight vectors: α ∈ R n i+1 Bias vector: γ ∈ R n i+1 Mean vector: where t bin i , t bn i , t lin i for i ∈ [d], t lin d+1 and t am d+1 are given in Table 1.
Intuitively, a LIN layer is a linear transformation. A BN layer following a LIN layer is used to standardize and normalize the output of the LIN layer. A BIN layer is used to binarize the real-numbered output vector of the BN layer. In this work, we consider the sign function which is widely used in BNNs to binarize real-numbered vectors. An ARGMAX layer follows a LIN layer and outputs the index of the largest entry as the predicted class which is represented by a one-hot vector. (In case there is more than one such entry, the first one is returned.) Formally, given a BNN N = (t 1 , · · · , t d , t d+1 ) and an input x ∈ B n 1 ±1 , N( x) ∈ B s is a one-hot vector in which the index of the non-zero entry is the predicated class.

Binary Decision Diagrams
A BDD [8] is a rooted acyclic directed graph where non-terminal nodes v are labeled by Boolean variables var(v) and terminal nodes (leaves) v are labeled with values val(v) ∈ B, referred to as the 1-leaf and the 0-leaf respectively. Each non-terminal node v has two outgoing edges: hi(v) meaning var(v) = 1 and lo(v) meaning var(v) = 0. We will also refer to hi(v) and lo(v) as the hi and lo children of v respectively. Moreover, assuming that x 1 , · · · , x m is the variable ordering, for each node v with var(v) = x i and each v ∈ {hi(v), lo(v)} with var(v ) = x j , we have i < j. In the graphical representation of BDDs, hi(v) and lo(v) are depicted by solid and dashed lines respectively. MTBDDs are a variant of BDDs in which the terminal nodes are not restricted to be 0 or 1. A BDD is reduced if it (1) has only one 1-leaf and one 0-leaf, (2) does not contain a node v such that hi(v) = lo(v), and (3) does not contain two distinct non-terminal nodes v and v such that var(v) = var(v ), hi(v) = hi(v ) and lo(v) = lo(v ). Hereafter, we assume that BDDs are reduced.
Bryant [8] showed that BDDs can serve as a canonical form of Boolean functions. Given a BDD over variables x 1 , · · · , x m , each non-terminal node v with var(v) = x i represents a Boolean function ). Operations on Boolean functions can usually be efficiently implemented via manipulating their BDD representations. A good variable ordering is crucial for the performance of BDD manipulations while the task of finding an optimal ordering for a function is NP-hard. To store and manipulate BDDs efficiently, nodes are stored in a hash table and recent computed results are stored in a cache to avoid duplicated computations. In this work, we will use some basic BDD operations such as ITE for If-Then-Else, Xor for exclusive-OR, Xnor for exclusive-NOR (i.e., a Xnor b = ¬(a Xor b)) and SatAll( f v ) for the set of all solutions of the Boolean formula f v . We denote by L(v) the set SatAll( f v ). For easy reference, more operations are given in Appendix A.1.

BDD4BNN Design
In this section, we first present an overview of our BDD-based quantitative analysis framework BDD4BNN, and then provide details of the key components.

BDD4BNN Overview
An overview of BDD4BNN is depicted in Figure 2. BDD4BNN comprises four main components: Region2BDD, BNN2CC, BDD Model Builder, and Query Engine. For a fixed BNN N = (t 1 , · · · , t d , t d+1 ) and a region R of the input space of N, BDD4BNN constructs the BDDs (G out i ) i∈[s] to encode the input-output relation of N in the region R, where the BDD G out i corresponds to the class i ∈ [s]. Technically, the region R is partitioned into s parts represented by the BDDs (G out i ) i∈ [s] . For each query of a property, BDD4BNN analyzes (G out i ) i∈[s] and outputs the query result. The general workflow of our approach is as follows. First, Region2BDD builds up a BDD G in R from the region R which represents the desired input space of N for analysis. Second, BNN2CC transforms each block of the BNN N into a set of cardinality constraints (CCs) similar to [42,6]. Third, BDD Model Builder builds the BDDs (G out i ) i∈[s] from all the cardinality constraints and the BDD G in R . Finally, Query Engine answers queries by analyzing the BDDs (G out i ) i∈[s] . Our Query Engine currently supports two types of application queries: robustness analysis and interpretability.
In the rest of this section, we first introduce the key sub-component CC2BDD, which provides encoding of cardinality constraints into BDDs. We then provide details of the components Region2BDD, BNN2CC, and BDD Model Builder. The Query Engine will be described in Section 4.

CC2BDD: Cardinality Constraints to BDDs
A cardinality constraint is a constraint of the form n j=1 j ≥ k over a vector x of Boolean variables with length n, where the literal j is either x j or ¬ x j for each j ∈ [n]. Note that constraints of the form n j=1 j > k, n j=1 j ≤ k and n j=1 j < k are equivalent to n j=1 j ≥ k + 1, n j=1 ¬ j ≥ n − k and n j=1 ¬ j ≥ n − k + 1, respectively. We assume that 1 (resp. 0) is a special cardinality constraint that always holds (resp. never holds).
To encode n j=1 j ≥ k as a BDD, we observe that all the possible solutions of n j=1 j ≥ k can be compactly represented by a BDD-like graph shown in Figure 3, where each node is labeled by a literal, and a solid (resp. dashed) edge from a node labeled by j means that the value of the literal j is 1 (resp. 0). Thus, each path from the 1 -node to the 1-leaf through the j -node (where 1 ≤ j ≤ n) captures a set of valuations where j followed by a (horizontal) dashed line is set to be 0 while j followed by a (vertical) solid line is set to be 1, and all the other literals which are not along the path can take arbitrary values. Clearly, for each of these valuations, there are at least k positive literals, hence the constraint n j=1 j ≥ k holds. Based on the above observation, we build the BDD for n j=1 j ≥ k using Algorithm 1. It builds a BDD for each node in Figure 3, row-by-row (the index i in Algorithm 1) and from right to left (the index j in Algorithm 1). For each node at the i-th row and j-th column, the label of the node must be the literal i+ j−1 . We build the BDD . Finally, we obtain the BDD G 1,1 that encodes the solutions of n j=1 j ≥ k. Lemma 1. For each cardinality constraint n j=1 j ≥ k, a BDD G with O((n − k) · k) nodes can be computed in O((n − k) · k) time such that L(G) is the set of all the solutions of n j=1 j ≥ k.

Region2BDD: Input Regions to BDDs
In this paper, we consider the following two types of input regions.
-Input region based on Hamming distance. For an input u ∈ B n 1 ±1 and an integer where HD( x, u) denotes the Hamming distance between x and u. Intuitively, R( u, r) includes the input vectors which differ from u by at most r positions.
-Input region with fixed indices. For an input u ∈ B n 1 ±1 and a set of indices I ⊆ [n 1 ], Intuitively, R( u, I) includes the input vectors which differ from u only at the indices from I.
Note that both R( u, n 1 ) and R( u, [n 1 ]) denote the entire input space B n 1 ±1 . Recall that each input sample is an element from B n 1 ±1 . To represent the region R by a BDD, we transform each value ±1 into a Boolean value 1/0. To this end, for each input Therefore, R( u, r) and R( u, I) will be represented by R( u (b) , r) and R( u (b) , I), respectively. The transformation functions t lin i , t bn i , t bin i and t am d+1 of the LIN, BN, BIN, and ARGMAX layers (cf. Table 1) will be handled accordingly. Note that for convenience, vectors over the Boolean domain B may be directly given by u or x when it is clear from the context. Region Encoding under Hamming distance. Given an input u ∈ B n 1 and an integer r, the region R( u, r) can be expressed by a cardinality constraint n 1 j=1 j ≤ r (which is equivalent to n 1 j=1 ¬ j ≥ n 1 − r), where for every j ∈ [n 1 ], j = x j if u j = 0, otherwise j = ¬ x j . For instance, consider u = (1, 1, 1, 0, 0) and r = 2, we have: Thus, R((1, 1, 1, 0, 0), 2) can be expressed by the cardinality constraint ¬ By Algorithm 1, the cardinality constraint of R( u, r) can be encoded by the BDD G in u,r , such that L(G in u,r ) = R( u, r). Following Lemma 1, we get that: Lemma 2. For every input region R given by an input u ∈ B n 1 and an integer r, a BDD G in u,r with O(r · (n 1 − r)) nodes can be computed in O(r · (n 1 − r)) time such that L(G in u,r ) = R( u, r). Region Encoding under fixed indices. Given an input u ∈ B n 1 and a set of indices

BNN2CC: BNNs to Cardinality Constraints
As mentioned before, to encode the BNN N = (t 1 , · · · , t d , t d+1 ) as BDDs, we transform the BNN N into cardinality constraints from which the desired BDDs (G out i ) i∈[s] are constructed. To this end, we first transform each internal block t i : B n i ±1 → B n i+1 ±1 into n i+1 cardinality constraints, each of which corresponds to one of the outputs of t i . Then we transform the output block t d+1 : B n d+1 ±1 → B s into s(s − 1) cardinality constraints, where one output class yields (s − 1) cardinality constraints.
For each vector-valued function t, we denote by t ↓ j the (scalar-valued) function returning the j-th entry of the output of t.
Transformation for internal blocks. Consider the internal block t i : , where 1 denotes the vector of 1's with the width n i .
Let C i, j be the following cardinality constraint: Transformation for the output block. For the output block t d+1 : . For every j ∈ [s] \ { j}, we define the cardinality constraint C d+1, j as follows: For each internal block t i : B n i ±1 → B n i+1 ±1 , we denote by BNN2CC(t i ) the cardinality constraints {C i,1 , · · · , C i,n i+1 }. For each output class j ∈ [s], we denote by BNN2CC j (t d+1 ) the cardinality constraints {C d+1,1 , · · · C d+1, j−1 , C d+1, j+1 , · · · , C d+1,s }. By applying the above transformation to all the blocks of the BNN N = (t 1 , · · · , t d , t d+1 ), we obtain its cardinality constraint form such that all the cardinality constraints in BNN2CC j (t d+1 ) hold under the valuation u. It is straightforward to verify:

BDD Model Builder
The construction of the BDDs (G out i ) i∈[s] from the BNN N (b) and the input region R is done iteratively throughout the blocks. Initially, the BDD for the first block is built, which can be seen as the input-output relation for the first internal block. In the i-th iteration, as the input-output relation of the first (i − 1) internal blocks has been encoded into the BDD, we compose this BDD with the BDD for the block t i which is built from its cardinality constraints t (b) i , resulting in the BDD for the first i internal blocks. Finally, we obtain the BDDs (G out i ) i∈[s] of the BNN N, with respect to the input region R. Design choice. There are several design choices for efficiency consideration which we discuss as follows. First of all, to encode the input-output relation of an internal block t i into BDD from its cardinality constraints t (b) i = {C i,1 , · · · , C i,n i+1 }, it amounts to compute And j∈[n i+1 ] CC2BDD(C i, j ). A simple and straightforward approach is to initially compute a BDD G = CC2BDD(C i,1 ) and then iteratively compute the conjunction G = And(G, CC2BDD(C i, j )) of G and CC2BDD(C i, j ) for 2 ≤ j ≤ n i+1 .
Alternatively, we use a divide-and-conquer strategy to recursively compute the BDDs for the first half and the second half of the cardinality constraints respectively, and then apply the AND-operation. Our preliminary experimental results show that the latter approach often performs better (about 2 times faster) than the former one, although they generate the same BDD.
Second, constructing the BDD directly from the cardinality constraints t (b) i = {C i,1 , · · · , C i,n i+1 } becomes prohibitively costly when n i and n i+1 are large, as the BDDs CC2BDD(C i, j ) for j ∈ [n i+1 ] need to consider all the inputs in B n i . To improve efficiency, we apply feasible input propagation. Namely, when we construct the BDD for the block t i+1 , we only consider its possible inputs with respect to the output of the block t i . Our preliminary experimental results show that the optimization could significantly improve the efficiency of the BDD construction.
Third, instead of encoding the input-output relation of the BNN N as a sole BDD or MTBDD, we opt to use a family of s BDDs (G out i ) i∈[s] , each of which corresponds to one output class of N. Recall that each output class i ∈ [s] is represented by (s − 1) cardinality constraints. Then, we can build a BDD G i for the output class i, similar to the BDD construction for internal blocks. By composing G i with the BDD of the entire internal blocks, we obtain the BDD G out i . Building a single BDD or MTBDD for the Algorithm 2: BDD Construction for BNNs , but our approach gives the flexibility especially when a specific target class is interested, which is common for robustness analysis.
Overall algorithm. The overall BDD construction procedure is shown in Algorithm 2. Given a BNN N = (t 1 , · · · , t d , t d+1 ) with s output classes and an input region R( u, τ), the algorithm outputs the BDDs (G out i ) i∈[s] , encoding the input-output relation of the BNN N with respect to the input region R( u, τ).
The procedure BNN2BDD first builds the BDD representation G in u,τ of the input region R( u, τ) and the cardinality constraints from BNN N (b) (Line 1). The first forloop builds a BDD encoding the input-output relation of the entire internal blocks w.r.t. G in u,τ . The second for-loop builds the BDDs (G out i ) i∈[s] , each of which encodes the input-output relation of the entire BNN for a class i ∈ [s] w.r.t. G in u,τ . The procedure Block2BDD receives the cardinality constraints {C m , · · · , C n }, a BDD G in representing the feasible inputs of the block and the block index i as inputs, and returns a BDD G. If i = d + 1, namely, the cardinality constraints {C m , · · · , C n } are from the output block, the resulting BDD G encodes the subset of G in u,τ that satisfy all the cardinality constraints {C m , · · · , C n }. If i d + 1, then the BDD G encodes the input-output relation of the Boolean function f m,n such that for every x i ∈ L(G in ), f m,n ( x i ) is the truth vector of the cardinality constraints {C m , · · · , C n } under the valuation x i . When m = 1 and n = n i+1 , f m,n is the same as Theorem 2. Given a BNN N with s output classes and an input region R( u, τ), we can compute s BDDs

Applications: Robustness Analysis and Interpretability
In this section, we present two applications within BDD4BNN, i.e., robustness analysis and interpretability of BNNs.

Robustness Analysis
Definition 2. Given a BNN N and an input region R( u, τ), the BNN is (locally) robust w.r.t. the region R( u, τ) if each sample x ∈ R( u, τ) is classified into the same class as the ground-truth class of u.
An adversarial example in the region R( u, τ) is a sample x ∈ R( u, τ) such that x is classified into a class, that differs from the ground-truth class of u.
As mentioned in Section 1, qualitative verification which checks whether a BNN is robust or not is insufficient in many practical applications. In this paper, we are interested in quantitative verification of robustness which asks how many adversarial examples are there in the input region of the BNN for each class. To answer this question, given a BNN N and an input region R( u, τ), we first obtain the BDDs (G out i ) i∈[s] by applying Algorithm 2 and then count the number of adversarial examples for each class in the input region R( u, τ). Note that counting adversarial examples amounts to computing |R( u, τ)| − |L(G out g )|, where g denotes the ground-truth class of u, and |L(G out g )| can be computed in time O(|G out g |). In some applications, more refined analysis is needed. For instance, it may be acceptable to misclassify a dog as a cat, but unacceptable to misclassify a tree as a car. This suggests that the robustness of BNNs may depend on the classes to which samples are misclassified. To capture this, we consider the notion of targeted robustness.
Definition 3. Given a BNN N, an input region R( u, τ), and the class t, the BNN is t-target-robust w.r.t. the region R( u, τ) if every sample x ∈ R( u, τ) is never classified into the class t. (Note that we assume that the ground-truth class of u differs from the class t.) The quantitative verification problem of t-target-robustness of a BNN asks how many adversarial examples in the input region R( u, τ) are misclassified to the class t by the BNN N. To answer this question, we first obtain the BDD G out t by applying Algorithm 2 and then count the number of adversarial examples by computing |L(G out t )|. Note that, if one wants to compute the (locally) maximal safe Hamming distance that satisfies a robustness property for an input sample (e.g., the proportion of adversarial examples is below a threshold), our framework can incrementally compute such a distance without constructing the BDD models of the entire BNN from scratch.  if Pr(R( u, r)) > , then Pr(R( u, r 1 )) ≤ and Pr(R( u, r )) > for r : r 1 < r < r; if Pr(R( u, r)) ≤ , then Pr(R( u, r 1 + 1)) > and Pr(R( u, r )) ≤ for r : r < r ≤ r 1 ; where Pr(R( u, r)) is the probability i∈[s].i g |L(G out i )| |R( u,r)| for g being the ground-truth class of u.
Algorithm 3 shows the procedure to incrementally compute the maximal safe Hamming distance for a given threshold ≥ 0, input region R( u, r), and ground-truth class g of u. Remark that Pr(R( u, r)) may not be monotonic w.r.t. the Hamming distance r.

Interpretability
In general, interpretability addresses the question of why some inputs in the input region are (mis)classified by the BNN into a specific class? We consider the interpretability of BNNs using two complementary explanations, i.e., prime implicant explanations and essential features.
Definition 5. Given a BNN N, an input region R( u, τ) and a class g, a prime implicant explanation (PI-explanation) of decisions made by the BNN N on the inputs L(G out g ) is a minimal set of literals { 1 , · · · , k } such that for every x ∈ R( u, τ), if x satisfies 1 ∧· · ·∧ k , then x is classified into the class g by the BNN N.
Intuitively, a PI-explanation { 1 , · · · , k } indicates that {var( 1 ), · · · , var( k )} are key features, namely, if fixed, the predication is guaranteed no matter how the remaining features change. Remark that there may be more than one PI-explanation for a set of inputs L(G out g ). When g is set to be the class of the benign input u, a PI-explanation on G out g suggests why these samples are classified into g by the BNN N.
Definition 6. Given a BNN N, an input region R( u, τ) and a class g, the essential features for the inputs L(G out g ) are literals { 1 , · · · , k } such that every x ∈ R( u, τ), if x is classified into the class g by the BNN N, then x satisfies 1 ∧ · · · ∧ k .
Intuitively, the essential features { 1 , · · · , k } denote the key features such that all samples x ∈ R( u, τ) that are classified into the class g by the BNN N must agree on these features. Essential features differ from PI-explanations, where the former can be seen as a necessary condition, while the latter can be seen as a sufficient condition.

Evaluation
We have implemented our framework as a prototype tool BDD4BNN based on the CUDD package [54]. BDD4BNN is implemented with Python as the front-end to preprocess BNNs and C++ as the back-end to perform the BDD encoding and analysis. In this section, we report the experimental results, including BDD encoding, robustness analysis, and interpretability. Because of space restriction, the results of BDD encoding and robustness analysis of BNNs with fixed indices are given in Appendix A.6.
Experimental Setup. The experiments were conducted on a machine with Intel Xeon Gold 5118 2.3GHz CPU, 64-bit Ubuntu 20.04 LTS operating systems, 128G RAM. Each BDD encoding executed on one core limited by 8-hour.
Benchmarks. We use the PyTorch (v1.0.1.post2) deep learning platform provided by NPAQ [6] to train and test BNNs. We trained 12 BNN models (P1-P12) with varying sizes using the MNIST dataset [34]. The MNIST dateset contains 70,000 gray-scale 28 × 28 images (60,000 for training and 10,000 for testing) of handwritten digits with 10 classes. In our experiments, we downscale the images (28 × 28) to some selected input size n 1 (i.e., the corresponding image is of the size √ n 1 × √ n 1 ) and then binarize the normalized pixels of the images.
Details of the BNN models are listed in Table 2, each of which has 10 classes (i.e., s = 10). Column 1 shows the name of the BNN model. Column 2 shows the architecture of the BNN model, where n 1 : · · · : n d+1 : s denotes that the BNN model has d + 1 blocks, n 1 inputs and s outputs; the i-th block for i ∈ [d + 1] has n i inputs and n i+1 outputs with n d+2 = s. Recall that each internal block has 3 layers while the output block has 2 layers. Therefore, the number of layers ranges from 5 to 14, the dimension of inputs ranges from 9 to 784, and the number of hidden neurons per linear layer ranges from 10 to 100. Column 3 shows the accuracy of the BNN model on the test set of the MNIST dataset. (We can observe that the accuracy increases with the size of inputs, the number of layers, and the number of hidden neurons per layer.) We randomly choose 10 images (shown in Figure 6 in Appendix) from the test set of the MNIST dataset (one image per class) to evaluate our approach.

Performance of BDD Encoding
We evaluate BDD4BNN on the BNNs listed in Table 2 using different input regions.
BDD encoding using full input space. We evaluate BDD4BNN on the BNNs (P1-P5), where B n 1 ±1 is used as the input region. The results are shown in Table 3, where |G| denotes the number of BDD nodes in the BDD manager. We can observe that both the execution time and the number of BDD nodes increase with the size of BNNs.
BDD encoding under Hamming distance. We evaluate BDD4BNN on the relatively large BNNs (P5-P12). In this case, an input region is given by one of the 10 images and a Hamming distance r ranging from 2 to 6. The average results are shown in Table 4, where [i] (resp. (i)) indicates the number of cases that BDD4BNN runs out of memory (resp. time). Overall, the execution time and the number of BDD nodes increase with r. BDD4BNN succeeded on all the cases when r ≤ 4, 75 cases out of 80 when r = 5, and 48 cases out of 80 when r = 6. We observe that the execution time and number of BDD nodes increase with the number of hidden neurons (P6 vs. P7, P8 vs. P9, and P11 vs. P12), while the effect of the number of layers is diverse (P6 vs. P8 vs. P10, and P7 vs. P9). From P9 and P10, we observe that the number of hidden neurons per layer is likely the key impact factor of the efficiency of BDD4BNN. Interestingly, our tool BDD4BNN works well on BNNs with large input sizes (i.e., on P11 and P12).
These results demonstrate the efficiency and scalability of BDD4BNN on BDD encoding of BNNs. We remark that, compared with the learning-based approach [50], our approach is considerably more efficient and scalable. For instance, the learning-based approach takes 403 seconds to encode a BNN with 64 input size, 5 hidden neurons, and 2 output size when r = 6, while ours takes about 3 seconds even for a larger network P5.

Robustness Analysis
We evaluate BDD4BNN on the robustness of BNNs, including robustness analysis under different input regions and maximal safe Hamming distance computing.
Robustness verification with Hamming distance. We evaluate BDD4BNN on BNNs (P7, P8, P9, and P11) using the 10 images. The input regions are given by the Hamming  distance r ranging from 2 to 4, resulting in 120 instances. To the best of our knowledge, NPAQ [6] is the only work that supports quantitative robustness verification of BNNs to which we compare BDD4BNN. Recall that NPAQ only provides PAC-style guarantees. Namely, it sets a tolerable error ε and a confidence parameter δ. The final estimated results of NPAQ have the bounded error ε with confidence of at least 1 − δ, i.e., In our experiments, we set ε = 0.8 and δ = 0.2, as done in [6]. The results on the average of the images are shown in Table 5. NPAQ ran out of time on 5 instances (which occur in P9 with r = 4 and P11 with r = 3 and r = 4), while BDD4BNN successfully verified all the 120 instances.   Table 4 and Table 5, we also found that most of the verification time is spent on BDD encoding while the rest is usually less than 10 seconds.
Details of robustness and targeted robustness. Figure 4(a) (resp. Figure 4(b) and Figure 4(c)) depicts the distributions of classes on P8 with Hamming distance r = 2 (resp. P8 with r = 3 and P11 with r = 2), where on the x-axis i = 0, · · · , 9 denotes the input region that is within the respective Hamming distance to the image of digit i (called i-region). We can observe that P8 is robust for the 0-region when r = 2 3 and robust for the 6-region when r = 2 and r = 3, but is not robust for the other regions. Most of the adversarial examples in the 1-region and 5-region are misclassified into the digit 3 by P8. P11 is not robust for the 1-region or the 5-region, but is robust for all the other regions. Though P8 and P11 are not robust on some input regions, indeed they are t-target-robust for many target classes t, e.g., P11 is t-target-robust for the 1-region when t 2, and the 5-region when t 3. (The raw data are given in Tables 8 and 9 in Appendix.) Quality validation of NPAQ.  ( = 0 and = 0.03). The initial Hamming distance r is 3. Intuitively, = 0 (resp. = 0.03) means that up to 0% (resp. 3%) samples in the input region can be adversarial. Table 6 shows the results, where columns SD and Time give the maximal safe Hamming distance and the execution time, respectively. BDD4BNN solved 74 out of 80 instances. (For the remaining 6 instances, BDD4BNN ran out of time or memory, but it was still able to compute a larger safe Hamming distance.) We can observe that the maximal safe Hamming distance increases with the threshold on several BNNs and input regions. We can also observe that P11 is more robust than others, which is consistent with their accuracies (cf. Table 2). Remark that SD = −1 indicates that the input image itself is misclassified.

Interpretability
To demonstrate the ability of BDD4BNN on interpretability, we consider the analysis of the BNN P12 and the image u of digit 1.
Essential features. For the input region given by the Hamming distance r = 4, we compute two sets of essential features for the inputs L(G out 2 ) and L(G out 5 ), i.e., the adversarial examples in the region R( u, 4) that are misclassified into the classes 2 and 5 respectively. The essential features are depicted in Figures 5(a) and 5(b), where black (resp. blue) color means that the value of the corresponding pixel is 1 (resp. 0), and yellow color means that the value of the corresponding pixel can take arbitrary values. Figure 5(a) (resp. Figure 5(b)) indicates that the inputs L(G out 2 ) (resp. L(G out 5 )) must agree on these black-and blue-colored pixels.
PI-explanations. For demonstration, we assume that the input region is given by the fixed set of indices I = {1, 2, · · · , 28} which denotes the first row of pixels of 28 × 28 images. We compute two PI-explanations of the inputs L(G out 2 ) and L(G out 5 ). The PIexplanations are depicted in Figures 5(c) and 5(d). Figure 5(c) (resp. Figure 5(d)) suggests that, by the definition of the PI-explanation, all the images in the region R( u, I) obtained by assigning arbitrary values to the yellow-colored pixels are always misclassified into the class 2 (resp. class 5), while changing one black-colored or blue-colored pixel would change the predication result since a PI-explanation is a minimal set of literals. In this section, we discuss the related work to BDD4BNN on qualitative/quantitative analysis and interpretability of DNNs. As there is a vast amount of literature regarding these topics, we will only discuss the most related ones.
Qualitative analysis of DNNs. For the verification of real-numbered DNNs, we broadly classify the existing approaches into three categories: (1) constraint solving based, (2) optimization-based, and (3) program analysis based.
The first class of approaches represents the early efforts which reduce to constraint solving. Pulina and Tacchella [47] verified whether the output of the DNN is within an interval by reducing to the satisfiability checking of a Boolean combination of linear arithmetic constraints via SMT solvers. Spurious adversarial examples can trigger refinements and retraining of DNNs. Katz et al. [29] and Ehlers [15] independently implemented two SMT solvers, Reluplex and Planet, for verifying properties of DNNs that are expressible with respective constraints. Recently, Reluplex was re-implemented in a new framework Marabou [30] with significant improvements.
For the second class of approaches which reduce to an optimization problem, Lomuscio and Maganti [37] verified whether some output is reachable from a given input region by reducing to mixed-integer linear programming (MILP) via optimization solvers. To speed up DNN verification via MILP solving, Cheng et al. [11] proposed heuristics for MILP encoding and parallelization of MILP solvers. Dutta et al. [13] proposed an algorithm to estimate the output region for a given input region. The algorithm iterates between a global search with MILP solving and a local search with gradient descent. Tjeng et al.
[57] proposed a tighter formulation for non-linearities in MILP and methods to improve performance. Recently, Bunel et al. [9] presented a branch and bound algorithm (BaB) to verify DNNs on properties expressible in Boolean formulas over linear inequalities. They claimed that both previous SAT/SMT and MILP-based approaches are its special cases. Convex optimization has also been used to verify DNNs with over-approximations [65, 14,67].
For the third class, researchers have adapted various methods from traditional static analysis to DNNs. A typical example is to use abstract interpretation, possibly aided with a refinement procedure to tighten approximations [18,52,53,35,2,51,36,68,58,59]. These methods vary in the abstract domain (e.g., box, zonotope, polytope, and star-set), efficiency, precision, and activation functions. (Remark that [52,53] considered floatingpoints instead of real numbers.) Another type is to compute convergent output bounds by exploring neural networks layer-by-layer. Huang et al. [25] proposed an exhaustive search algorithm with an SMT-based refinement. Later, the search problem was solved via Monte-Carlo tree search [64,66]. Xiang et al. [63] proposed to approximate the bounds based on the linear approximations for the neurons and Lipschitz constants [22]. Wang et al. [60] presented symbolic interval analysis to tighten approximations. Recently, two abstraction-based frameworks have been proposed [3,16] which aim to reduce the size of DNNs, making them more amenable to verification.
Quantitative analysis of DNNs. Comparing to the qualitative analysis, the quantitative analysis of neural networks is currently very limited. Two sampling-based approaches were proposed to certify the robustness of adversarial examples [61,5], which require only blackbox access to the models, hence can be applied on both DNNs and BNNs. Yang et al.
[68] proposed a spurious region-guided refinement approach for real-numbered DNN verification. The quantitative robustness verification is achieved by over-approximating the Lebesgue measure of the spurious regions. The authors claimed that it is the first work to quantitative robustness verification of DNNs with soundness guarantee.
Following the SAT-based qualitative analysis of BNNs [42,41], SAT-based quantitative analysis approaches were proposed [6,43,19] for verifying robustness and fairness, and assessing heuristic-based explanations of BNNs. In particular, approximate SAT model-counting solvers are utilized. As demonstrated in Section 5, our BDD-based approach is considerably more accurate and efficient than the SAT-based one [6]. In general, we remark that the BDD construction is computationally expensive but the follow-up analysis is often much more efficient, while the SAT encoding is efficient (polynomial-time) but SAT queries are often computationally expensive ( P-hard). The computational cost of our approach is more dependent on the number of neurons per linear layer but less on the number of layers, while the computational cost of the SAT-based approach [6] is dependent on both of them.
Shih et al.
[50] proposed a BDD-based approach to tackle BNNs, similar to our work, in spirit. In this BDD learning-based approach, membership queries are implemented by querying the BDD for each input, equivalence queries are implemented by transforming the BDD and BNN to two Boolean formulas, and checking the equivalence of two Boolean formulas under the input region (in a Boolean formula) via SAT solving. This construction requires n equivalence queries and 6n 2 + n · log(m) membership queries, where n (resp. m) is the number of nodes (resp. variables) in the final BDD. Compared with this approach, our approach is able to handle much larger BNNs than theirs.
Interpretability of DNNs. Though interpretability of DNNs is crucial for explaining predictions, it is very challenging to tackle due to the blackbox nature of DNNs. There is a large body of work on the interpretability of DNNs (cf. [24,39] for a survey). Almost all the existing approaches are heuristic-based and restricted to finding explanations that are local in an input region. Some of them tackle the interpretability of DNNs by learning an interpretable model, such as binary decision trees [17,69] and finite-state automata [62]. In contrast to ours, they target at DNNs and only approximate the original model in the input region. The BDD-based approach [50] mentioned above has been used to compute PI-explanation while essential features were not considered therein.

Conclusion
In this paper, we have proposed a novel BDD-based framework for the quantitative verification of BNNs. We implemented our framework as a prototype tool BDD4BNN and conducted extensive experiments on 12 BNN models with varying sizes and input regions. Experimental results demonstrated that BDD4BNN is more scalable than the existing BDD-learning based approach, and significantly efficient and accurate than the existing SAT-based approach NPAQ. This work represents the first, but a key step of the long-term program to develop an efficient and scalable BDD-based quantitative analysis framework for BNNs.    Table 7 provides the BDD operations used in this work. We denote by op(v, v ) the operation Apply(v, v , op).

A.2 Proof of Proposition 1
Consider the internal block t i : where W k, j , b j , µ j , α j , σ j , γ j for k ∈ [n i ] are constants. Therefore, for every x ∈ B n i , we have: otherwise.

A.4 Explanation of Algorithm 2
Given a BNN N = (t 1 , · · · , t d , t d+1 ) with s output classes and an input region R( u, τ), Algorithm 2 outputs the BDDs (G out i ) i∈[s] , encoding the input-output relation of the BNN N w.r.t. the input region R( u, τ).
In detail, it first builds the BDD representation G in u,τ of the input region R( u, τ) and the cardinality constraints from BNN N (b) (Line 1). The first for-loop. It builds a BDD encoding the input-output relation of the entire internal blocks w.r.t. G in u,τ . It first invokes the procedure Block2BDD(t b i , G in , i) to build a BDD G encoding the input-output relation of the i-th block t b i w.r.t. L(G in ) (Line 4). G in is the set of feasible inputs of the block t b i , which is also the set of feasible outputs of the (i − 1)-th block t b i−1 (the input region G in u,τ when i = 1). By doing so, we have: From the BDD G , we compute the feasible outputs G in of the block t (b) i by existentially quantifying all the input variables x i of the block t (b) i (Line 5). The BDD G in serves as the set of feasible inputs of the block t (b) i+1 at the next iteration. We next assign G to G if the current block is the first internal block (i.e., i = 1), otherwise we compute the relational product of G and G, the resulting BDD G encodes the input-output relation of the first i internal blocks w.r.t. G in u,τ (Line 6), namely, At the end of the first for-loop, we obtain the BDD G encoding the input-output relation of the entire internal blocks and its feasible outputs G in w.r.t. G in u,τ , namely, The second for-loop. It builds the BDDs (G out i ) i∈[s] , each of which encodes the inputoutput relation of the entire BNN and a class i ∈ [s] w.r.t. G in u,τ . For each i ∈ [s], it first Fig. 6: 10 randomly chosen images used to evaluate our approach builds a BDD G i encoding the input-output relation of the output block for the class i by invoking Block2BDD(t (b) d+1↓i , G in , d + 1) (Line 8). By computing the relational product of the BDDs G i and G, we obtain the BDD G out i . Recall that the BDD G encodes the input-output relation of the entire internal blocks w.r.t. the input region G in u,τ . Thus, an input x ∈ R( u, τ) is classified into the class i by the BNN N iff x (b) ∈ L(G out i ). The procedure Block2BDD. It receives the cardinality constraints {C m , · · · , C n } (note that indices matter), a BDD G in encoding feasible inputs and the block index i as inputs, and returns a BDD G.
-If i = d + 1, namely, the cardinality constraints {C m , · · · , C n } are from the output block, the resulting BDD G encodes the subset of G in u,τ that satisfy all the cardinality constraints {C m , · · · , C n }.
-If i d+1, then the BDD G encodes the input-output relation of the Boolean function f m,n such that for every x i ∈ L(G in ), f m,n ( x i ) is the truth vector of the cardinality constraints {C m , · · · , C n } under the valuation x i . When m = 1 and n = n i+1 , f m,n is the same as In detail, the procedure Block2BDD computes the desired BDD in a binary search fashion.
-If m = n, it first builds the BDD G 1 for the cardinality constraint C n (Line 13) such that L(G 1 ) is the set of solutions of C n . Then it computes the conjunction G of the BDDs G 1 and G in . Recall that G in is the set of feasible input of the i-th block t (b) i . Thus, L(G) ⊆ L(G in ) is the set of feasible inputs that satisfy C n . If i d + 1, the BDD G is transformed into the BDD Xnor( x i+1 m , G). This step encodes the constraint x i+1 m ⇔ C n in the BDD Xnor( x i+1 m , G), namely, for every ( u, b) ∈ L(Xnor( x i+1 m , G)), b is the truth of C n under the valuation u. Remark that x i+1 m is a new Boolean variable introduced into the BDD Xnor( x i+1 m , G). -If m n, we recursively build the BDDs G 1 and G 2 for {C m , · · · , C n−m 2 +m } and {C n−m 2 +m+1 , · · · , C n } and compute the conjunction of G 1 and G 2 . Thus, if i = d + 1, L(G) ⊆ L(G in ) is the set of all the feasible inputs that satisfy the cardinality constraints {C m , · · · , C n }. If i d + 1, for every ( u, x) ∈ L(G), x is the truth vector of the constraints {C m , · · · , C n } under the valuation u.
A.5 Distributions of Classes on P8 and P11 with Hamming Distance r = 2, 3, 4 Table 8 shows the number of samples for each input region and class on P8, while Table 9 shows the number of samples for each input region and class on P11, where the input regions are given by Hamming Distance r = 2, 3, 4.