Model-based feature selection for neural networks: A mixed-integer
  programming approach
In this work, we develop a novel input feature selection framework for
ReLU-based deep neural networks (DNNs), which builds upon a mixed-integer
optimization approach. While the method is generally applicable to various
classification tasks, we focus on finding input features for image
classification for clarity of presentation. The idea is to use a trained DNN,
or an ensemble of trained DNNs, to identify the salient input features. The
input feature selection is formulated as a sequence of mixed-integer linear
programming (MILP) problems that find sets of sparse inputs that maximize the
classification confidence of each category. These ''inverse'' problems are
regularized by the number of inputs selected for each category and by
distribution constraints. Numerical results on the well-known MNIST and
FashionMNIST datasets show that the proposed input feature selection allows us
to drastically reduce the size of the input to $\sim$15\% while maintaining a
good classification accuracy. This allows us to design DNNs with significantly
fewer connections, reducing computational effort and producing DNNs that are
more robust towards adversarial attacks.
Shudian Zhao, Calvin Tsay, Jan Kronqvist
# Model-based feature selection for neural networks: A mixed-integer programming approach+

###### Abstract

In this work, we develop a novel input feature selection framework for ReLU-based deep neural networks (DNNs), which builds upon a mixed-integer optimization approach. While the method is generally applicable to various classification tasks, we focus on finding input features for image classification for clarity of presentation. The idea is to use a trained DNN, or an ensemble of trained DNNs, to identify the salient input features. The input feature selection is formulated as a sequence of mixed-integer linear programming (MILP) problems that find sets of sparse inputs that maximize the classification confidence of each category. These "inverse" problems are regularized by the number of inputs selected for each category and by distribution constraints. Numerical results on the well-known MNIST and FashionMNIST datasets show that the proposed input feature selection allows us to drastically reduce the size of the input to \(\sim\)15% while maintaining a good classification accuracy. This allows us to design DNNs with significantly fewer connections, reducing computational effort and producing DNNs that are more robust towards adversarial attacks.

_Keywords--_ Mixed-integer programming, Deep neural networks, Feature selection, Sparse DNNs, Model reduction.

## 1 Introduction

Over the years, there has been an active interest in algorithms for training sparse deep neural networks (DNNs) or sparsifying trained DNNs. By sparsifying a DNN we mean removing some connections (parameters) in the network, which can be done by setting the corresponding weights to zero. Examples of algorithms for sparsifying or training sparse networks include, dropout methods [8; 14; 15; 31], optimal/combinatorial brain surgeon [13; 36], optimal brain damage [19], and regularization based methods [21; 23; 33]. Carefully sparsifying the network, _i.e.,_ reducing the number of parameters wisely, has shown to reduce
over-fitting and improve overall generalizability [14, 18]. This paper focuses on feature selection for DNNs, which can also be interpreted as "sparsifying" the first/input layer, and we show that we can significantly reduce the number of parameters, _i.e._, non-zero weights, while keeping a good accuracy. Throughout the paper we focus on image classification, but the framework is general.

This work focuses on finding the salient input features for classification using a DNN. We hypothesize that the number of inputs to DNNs for classification can often be greatly reduced by a "smart" choice of input features while keeping a good accuracy (_i.e.,_ feature selection). We build the hypothesis on the assumption that not all inputs, or pixels, will be equally important. Reducing the number of inputs/parameters has the potential to: i) reduce over-fitting, ii) give more robust DNNs that are less sensitive for adversarial attacks (fewer degrees of freedom for the attacker), and iii) reduce computational complexity both in training and evaluating the resulting classifier (fewer weights to determine and fewer computational operations to evaluate the outputs). The first two are classical focus areas within artificial intelligence (AI), and the third is becoming more important with an increasing interest in so-called green AI [26]. Most strategies for feature selection can be grouped as either _filter_ methods, which examine the data, _e.g.,_ for correlation, and _wrapper_ methods, which amount to a guided search over candidate models [4, 20]. Feature selection can also be incorporated directly into model training using _embedded_ methods, _e.g._, regularization. This paper and the numerical results are intended as a proof of concept to demonstrate mixed-integer linear programming (MILP) as an alternative technology for extracting the importance of input features from DNNs. Input feature selection is an active research area, _e.g.,_ see the review papers [9, 37], and a detailed comparison to state-of-the-art methods is not within the scope of this paper.

Our proposed method leverages trained models that achieve desirable performance, and attempts to select a feature set that replicates the performance using mixed-integer programming. We build on the idea that, given a relatively well-trained DNN, we can analyze the DNN in an inverse fashion to derive information about the inputs. Specifically, to determine the most important inputs, or pixels in the case of image classification, for a given label, we solve an optimization problem that maximizes the classification confidence of the label with a cardinality constraint on the number of non-zero inputs to the DNN. We consider this as an "inverse problem," as the goal is to determine the DNN inputs from the output. We additionally propose some input distribution constraints to make the input feature selection less sensitive to errors in the input-output mapping of the DNN. We only consider DNNs with the rectified linear unit (ReLU) activation function, as it enables the input feature selection problem to be formulated as a MILP problem [7, 22]. However, the framework can be easily generalized to CNN architectures and other MIP representable activation functions, e.g., max pooling and leaky ReLU

Optimizing over trained ReLU-based DNNs has been an active research topic in recent years, and has a wide variety of applications including verification [2, 22, 28], lossless compression [27], and surrogate model optimization [11, 35]. There even exists software, such as OMLT [3], for directly incorporating ReLU DNNs into general optimization models. Optimizing over a trained ReLU-based DNN through the MILP encoding is not a trivial task, but significant progress has been made in terms of strong formulations [1, 16, 29], solution methods [5, 25], and techniques for deriving strong valid inequalities [1, 2]. In combination with the remarkable performance of state-of-the-art MILP solvers, optimization over DNNs appears computationally tractable (at least for moderate size DNNs). This work builds upon recent optimization advancements, as reliably optimizing over ReLU-based DNNs is
a key component in the proposed method.

The paper is structured as follows. Section 2 first describes the MILP problem to determine which inputs maximize the classification confidence. Some enhancements for the input selection are presented, and the complete input selection algorithm is presented in Section 2.4. Numerical results are presented in Section 3, where we show that we can obtain a good accuracy when downsizing the input to 15% by the proposed algorithm, and that the resulting DNNs are more robust towards adversarial attacks in the \(\ell_{\infty}\) sense. Section 4 provides some conclusions.

## 2 Input feature selection algorithm

Our feature selection strategy is based on the idea of determining a small optimal subset of inputs, or pixels for the case of image classification, that are allowed to take non-zero values to maximize the classification confidence for each label using a pre-trained DNN. By combining the optimal subsets for each label, we can determine a set of salient input features. These input features can be considered as the most important for the given DNN, but we note that the DNN might not offer a perfect input-output mapping. To mitigate the impact of model errors, we propose a technique of using input distribution constraints to ensure that the selected input features are to some extent distributed over the input space. This framework could easily be extended to use optimization over an ensemble DNN model [32] for input selection, where the inputs would be selected such that the ensemble classification confidence is maximized. While using DNN ensembles can further mitigate the effect of errors in individual DNNs, our initial tests did not indicate clear advantages of using DNN ensembles for this purpose.

Here we focus on fully connected DNNs that classify grayscale images into 10 categories. While this setting is limited, the proposed method is applicable to classification of RGB images, and other classification problems in general. The input features are scaled between 0 and 1, with 255 (white) corresponding to 1 and 0 remaining black. We start by briefly reviewing MILP encoding of DNNs in the next subsection, and continue with more details on the input feature selection in the following subsections.

### Encoding DNNs as MILPs

In a fully connected ReLU-based neural network, the \(l\)-th layer with input \(x^{l}\) and output \(x^{l+1}\) is described as

\[x^{l+1}=\max\{0,W^{l}x^{l}+b^{l}\},\]

where \(W^{l}\in\mathbb{R}^{n_{l+1}\times n_{l}}\) is the weight matrix and \(b^{l}\in\mathbb{R}^{n_{l+1}}\) is the bias vector.

The input-output mapping of the ReLU activation function is given by a piece-wise linear function, and is mixed-integer representable [30]. There are different formulations for encoding the ReLU activation function using MILP, where the big-M formulation [7, 22] was the first presented MILP encoding and remains a common approach. For the \(i\)-th ReLU node at a fully connected layer with input \(x^{l}\), the big-M formulation for input-output
relation is given by

\[\begin{gathered}(w_{i}^{l})^{\top}x^{l}+b_{i}^{l}\leq x_{i}^{l+1},\\ (w_{i}^{l})^{\top}x^{l}+b_{i}^{l}-(1-\sigma)LB_{i}^{l+1}\geq x_{i}^{l +1},\\ x_{i}^{l+1}\leq\sigma UB_{i}^{l+1},\\ \sigma\in\{0,1\},\ x_{i}^{l+1}\geq 0,\end{gathered} \tag{1}\]

where \(w_{i}^{l}\) is the \(i\)-th row vector of \(W^{l}\), \(b_{i}^{l}\) is the \(i\)-th entry of \(b^{l}\), \(LB_{i}^{l+1}\) and \(UB_{i}^{l+1}\) are upper and lower bounds on the pre-activation function over \(x_{i}^{l+1}\), such that \(LB_{i}^{l+1}\leq(w_{i}^{l})^{\top}x^{l}+b_{i}^{l}\leq UB_{i}^{l+1}\).

The big-M formulation is elegant in its simplicity, but it is known to have a weak continuous relaxation which may require the exploration of a huge number of branch-and-bound nodes in order to solve the problem [1]. Anderson et al. [1] presented a so-called extended convex hull formulation, which gives the strongest valid convex relaxation of each individual node, and a non-extended convex hull formulation. Even though the convex hull is the strongest formulation for each individual node, it does not in general give the convex hull of the full input-output mapping of the DNN. Furthermore, the convex hull formulation results in a large problem formulation that can be computationally difficult to work with. The class of partition-based, or \(P\)-split, formulations, was proposed as an alternative formulation with a stronger continuous relaxation than big-M and computationally cheaper than the convex hull [16, 29]. Computational results in [17, 29] show that the partition-based formulation often gives significant speed-ups compared to the big-M or convex hull formulations. Here, we do not focus on the computational efficiency of optimizing over ReLU-based DNNs, and, for the sake of clarity, we use the simpler big-M formulation (1). In fact, for the problems considered in this work, the computational time to solve the MILP problems did not represent a limiting factor (big-M tends to actually perform relatively well for simple optimization problems). But, alternative/stronger formulations could directly be used within our framework.

### The optimal sparse input features (OSIF) problem

With the MILP encoding of the DNN, we can rigorously analyze the DNN and find extreme points in the input-output mapping. Recently Kronqvist et al. [17] illustrated a similar optimal sparse input features (OSIF) problem. This problem aims to maximize the probability of at most \(\aleph\) non-zero input features being classified with a certain label \(i\) for a given trained DNN. The problem is formed by encoding the ReLU activation function for each hidden node by MILP. Instead of the softmax function, the objective function is \(x_{i}^{L}\), where \(x^{L}\) is the output vector, thus maximizing the classification confidence of label \(i\).
Using the big-M formulation, the OSIF problem can be stated as

max \[x_{i}^{L}\] (2a) s.t. \[W^{l}x^{l}+b^{l}\leq x^{l+1},\ \forall l\in[L-1],\] (2b) \[W^{l}x^{l}+b^{l}-\text{diag}(LB^{l+1})(\mathbf{1}-\sigma^{l+1}) \geq x^{l+1},\ \forall\ l\in[L-1],\] (2c) \[x^{l}\leq\text{diag}(UB^{l})\sigma^{l},\ \sigma^{l}\in\{0,1\}^{n_{l}}, \forall\ l\in\{2,\ldots,L\},\] (2d) \[x^{L}=W^{L-1}x^{L-1}+b^{L-1},x^{L}\in\mathbb{R}^{10},\] (2e) \[x^{l}\in\mathbb{R}_{+}^{n_{l}},\forall\ l\in[L-1],\] (2f) \[y\geq x^{1},y\in\{0,1\}^{n_{1}},\] (2g) \[\mathbf{1}^{\top}y\leq\aleph_{0},\] (2h)

where \(n_{1}\) is the size of the input data, \(x^{L}\in\mathbb{R}^{10}\) is the output, \(L\) is the number of layers, \(LB^{l}\) and \(UB^{l}\) are the bounds on \(x^{l}\), \(\mathbf{1}\) denotes the all-ones vector, \(\text{diag}(\cdot)\) denote the matrix with \(\cdot\) on the diagonal and \(0\) on other entries. Eq. (2g) and (2h) describe the cardinality constraint \(\|x^{1}\|_{0}\leq\aleph_{0}\), which limits the number of selected inputs. Fig 1 shows some example results of solving problem (2) with \(\aleph\in\{10,20\}\) and class \(i\in\{0,1\}\). Given a larger cardinality number \(\aleph\), the latter is more visually recognizable from the selected pixels.

### Input distribution constraints

To extract information across the whole input image, we propose to add constraints to force the selected pixels to be distributed evenly across some pre-defined partitioning of the input space. Forcing the selected pixels to be more spread-out may also mitigate the effect of inaccuracy of the DNN used in the OSIF problem, _e.g.,_ by preventing a small area to be given overly high priority. There are various ways to partition the input variables. In this

Figure 1: The results of OSIF for class \(0\) and \(1\) on MNIST (\(\aleph_{0}\in\{10,20\}\)). Note that the selected pixels are white.
paper, we focus on image classification problems with square images as input. Furthermore, we assume that the images are roughly centered. Therefore, we denote each input variable as matrix \(X\in\mathbb{R}^{n\times n}\) and define the partition as \(k^{2}\) submatrices of equal size, _i.e.,_\(X^{ij}\in\mathbb{R}^{\frac{n}{k}\times\frac{n}{k}}\) for \(i,j\in[k]\). For instance, given \(n\) is even and \(k=2\), a natural partition of the matrix is

\[X=\begin{pmatrix}X^{11}&X^{12}\\ X^{21}&X^{22}\end{pmatrix}.\]

In this way, we denote \(x^{1}:=\text{vec}(X)\) the input data and \(n_{1}:=n^{2}\) the size of the input data, then \(I_{ij}\) is the index set for entries mapped from \(X_{ij}\)

\[I_{ij}=\{(i_{1}-1)n+i_{2}\mid i_{1}\in\{(i-1)\frac{n}{k}+1,\ldots,i\frac{n}{k} -1\},i_{2}\in\{(j-1)\frac{n}{k}+1,\ldots,j\frac{n}{k}-1\}\}.\]

We denote the collection of index sets for the partition as \(\mathcal{I}:=\{I_{i,j}\}_{\forall i,j\in[k]}\). To limit the number of pixels selected from each box for each category we add the following constraints

\[\lfloor\frac{\aleph_{0}}{k^{2}}\rfloor\leq\sum_{i\in I_{t}}y_{i}\leq\lceil \frac{\aleph_{0}}{k^{2}}\rceil,\forall I_{t}\in\mathcal{I}, \tag{3}\]

The constraint (3) forces the pixels to spread evenly between all partitions, while allowing some to contain one more selected pixel for each category.

To illustrate on the impact on the distribution constraints, Fig. 2 compares the selected pixels for MNIST with \(k\in\{1,2\}\). Compared to the result without distribution constraints (equivalent to \(k=1\)), pixels selected with \(k=2\) are more scattered over the whole image and are more likely to identify generalizable distinguishing features of the full input data, assuming the dataset has been pre-processed for unused areas of the images matrices.

### Controlling the number of selected features

Repeatedly solving the OSIF problem (2) for each label, _i.e.,_\(i\in\{0,\ldots,9\}\), and taking the union of all selected pixels does not give us full control of the total number of selected pixels. Specifically, some of the pixels can be selected by the OSIF problem for multiple classes, resulting in fewer combined pixels (the union of selected subsets) than an initial target.

Therefore, we present an approach to control the number of selected inputs, which we use in the proposed MILP-based feature selection algorithm. The main idea is to allow freedom

Figure 2: Optimal input features for MNIST (\(\aleph=50\))
over features already selected by previous classes in the current OSIF problem and adjust the cardinality constraint (2h) to

\[\sum_{i\in[n_{1}]\setminus J}y_{i}\leq\aleph_{0}, \tag{4}\]

where \(J\) is the index set for input features selected by previous models. Similarly, constraints (3) are adjusted as

\[\lfloor\frac{\aleph_{0}}{k^{2}}\rfloor\leq\sum_{i\in I_{t}\setminus J}y_{i}\leq \lceil\frac{\aleph_{0}}{k^{2}}\rceil,\forall I_{t}\in\mathcal{I}. \tag{5}\]

Finally, we formulate the OSIF problem with input distribution and total number control as

\[OSIF(\mathcal{M},i,\aleph_{0},\mathcal{I},J)=\text{argmax }\{x_{i}^{L}\ \mid(\ref{eq:2b})-(\ref{eq:2g}),(\ref{eq:4}),(\ref{eq:5})\}. \tag{6}\]

Based on the described techniques, we introduce the input feature selection algorithm, which is presented as pseudo code in Alg. 1.

```
Data: the number of features \(\aleph\), a trained DNN \(\mathcal{M}\), a matrix partition set \(\mathcal{I}\), class set \(\mathcal{C}\); Input:\(J\leftarrow\emptyset\); Output: Index set \(J\); \(\aleph_{0}\leftarrow\aleph/10\); for\(i\in\mathcal{C}\)do \(x\gets OSIF(\mathcal{M},i,\aleph_{0},\mathcal{I},J)\) ; # Eq. (7) \(J\gets J\cup\{s\mid x_{s}^{1}=0,s\in[n_{1}]\}\);  end for
```

**Algorithm 1**MILP-based feature selection (MILP-based selection)

## 3 Computational results

In this paper, we focus on image classification problems for the MNIST [6] and FashionMNIST [34] datasets. Both datasets consist of a training set of 60,000 examples and a test set of 10,000 examples. Each sample image in both datasets is a \(28\times 28\) grayscale image associated with labels from 10 classes. MNIST is the dataset of handwritten single digits between 0 and 9. FashionMNIST is a dataset of Zalando's article images with 10 categories of clothing. There is one fundamental difference between the two data sets, besides FashionMNIST being a somewhat more challenging data set for classification. In MNIST there are significantly more pixels that do not change in any of the training images, or only change in a few images, compared to FashionMNIST. The presence of such "dead" pixels is an important consideration for input feature selection.

Image preprocessing and training DNNs are implemented in PyTorch [24], and the MILP problems are modeled and solved by Gurobi through the Python API [12]. We trained each DNN with 2 hidden layers of the same size.
### Accuracy of DNNs with sparse input features

The goal is to illustrate that Alg. 1 can successfully identify low-dimensional salient input features. We chose to focus on small DNNs, as DNN \(2\times 20\) can already achieve an accuracy of 95.7% (resp. 86.3%) for MNIST (resp. FashionMNIST) and larger DNNs did not give clear improvements for the input selection. For such models, the computational cost of solving the MILPs is low1.

Footnote 1: On a laptop with a 10-core CPU, Gurobi can solve instances with \(\aleph=100\) and a DNN of \(2\times 20\) under 15 seconds. However, previous research [1, 29] has shown that significant speed-ups can be obtained by using a more advanced MILP approach.

Table 1 and Table 2 present the accuracies of DNNs with sparse input features on MNIST and FashionMNIST. It is possible to obtain a much higher accuracy by considering a more moderate input reduction (about 0.5% accuracy drop with 200 -300 inputs), but this defeats the idea of finding low dimensional salient features. For grayscale input image of \(28\times 28\), we select at most 15% input features and present the results with \(\aleph\in\{50,100\}\).

#### 3.1.1 MILP-based feature selection

Table 1 compares the accuracy of classification models with different architectures, _i.e.,_ with \(2\times 20\) vs. with \(2\times 40\). We select sparse input features by Alg. 1 with OSIF models \(2\times 10\) and \(2\times 20\). Since the distribution constraints are supposed to force at least one pixel selected in each submatrix, we select partition number \(k\in\{1,2\}\) and \(k\in\{1,2,3\}\) for instances with \(\aleph=50\) and \(\aleph=100\) respectively.

First, we investigate the effect of the distribution constraints. Table 1 and Table 2 both show that the accuracy increases when adding the distribution constraints (_i.e._, from \(k=1\) to \(k=2\)) for \(\aleph\in\{50,100\}\). However, the distribution constraints become less important as the number of selected features \(\aleph\) increases; the best \(k\) for MNIST and FashionMNIST varies for \(\aleph=100\) (noting that the choice of \(k\) also affects accuracy less). For MNIST with \(\aleph=100\), the accuracy of instances drops slightly as the \(k\) increases from 2 to 3, while

\begin{table}
\begin{tabular}{c|c c c|c c c} \multicolumn{8}{c}{DNNs of \(2\times 20\) hidden layers} \\ \hline \(\aleph\) & OSIF Model & \(k\) & Acc. & OSIF Model & \(k\) & Acc. \\ \hline \multirow{2}{*}{50} & \(2\times 10\) & 1 & 80.6\% & \(2\times 20\) & 1 & 80.5\% \\  & \(2\times 10\) & 2 & **85.3\%** & \(2\times 20\) & 2 & **86.6\%** \\ \hline \multirow{2}{*}{100} & \(2\times 10\) & 1 & 88.8\% & \(2\times 20\) & 1 & 89.0\% \\  & \(2\times 10\) & 2 & **91.2\%** & \(2\times 20\) & 2 & **90.6\%** \\  & \(2\times 10\) & 3 & 89.3\% & \(2\times 20\) & 3 & 89.3\% \\ \hline \multicolumn{8}{c}{DNNs of \(2\times 40\) hidden layers} \\ \hline \(\aleph\) & OSIF Model & \(k\) & Acc. & OSIF Model & \(k\) & Acc. \\ \hline \multirow{2}{*}{50} & \(2\times 10\) & 1 & 83.1\% & \(2\times 20\) & 1 & 82.7\% \\  & \(2\times 10\) & 2 & **87.6\%** & \(2\times 20\) & 2 & **89.2\%** \\ \hline \multirow{2}{*}{100} & \(2\times 10\) & 1 & 91.4\% & \(2\times 20\) & 1 & 91.4\% \\  & \(2\times 10\) & 2 & **93.4\%** & \(2\times 20\) & 2 & **92.9\%** \\ \multicolumn{8}{c}{} \\ \(2\times 10\) & 3 & 92.3\% & \(2\times 20\) & 3 & 91.8\% \\ \hline \end{tabular}
\end{table}
Table 1: Accuracy of DNNs of different architectures with sparse input features selected by Alg. 1 on MNIST
using input features selected with \(k=3\) leads to slightly higher accuracy for FashionMNIST with \(\aleph=100\). One reason behind this difference could be that input pixels of MNIST and FashionMNIST are activated in different patterns. In MNIST, there are more active pixels in the center, and the peripheral pixels stay inactive over the training data. Given the distribution constraints with \(k=2\) (see Fig. 2(b)), the selected pixels stay away from the peripheral area. However, when \(k=3\) (see Fig. 2(c)), more pixels are forced to be chosen from the upper right-hand and bottom left-hand corners. In contrast, the active pixels are more evenly spread across the full input in FashionMNIST. Hence, as \(k\) increases (see Fig. 2(e) and Fig. 2(f)), the active pixels remain well-covered by the evenly scattered selected pixels.

Next, Table 1 and Table 2 also compare OSIF using different DNN architectures, _i.e._, \(2\times 10\) and \(2\times 20\). The accuracy is 94% (resp. 84%) for the former and 96% (resp. 86%) for the latter for MNIST (resp. FashionMNIST). The results show that even using a simple DNN for feature selection using OSIF can produce feature sets that achieve good accuracy, when appropriately large classification models are trained on the selected features. For MNIST (see Table 1), the former model has accuracy at most 2 points worse than the latter model when \(\aleph=50\). When \(\aleph=100\), both models achieve similar levels of performance. As for Fashion, the difference is at most 1 point for \(\aleph\in\{50,100\}\). Hence, we cannot observe a clear difference between the two OSIF DNN models in terms of feature selection quality.

Finally, we would like to make a brief remark on the improvement in accuracy by increasing the size of the DNNs for classification, _i.e.,_ from \(2\times 20\) to \(2\times 40\), for both MNIST and FashionMNIST. Unsurprisingly, using larger DNNs results in overall higher accuracy. More importantly, the performance of the proposed input feature selection seems to be robust toward the final architecture. For both architectures, we observe a similarly slight reduction in classification accuracy related to the reduction in number of input features (pixels).

\begin{table}
\begin{tabular}{c|c c c|c c c} \multicolumn{8}{c}{DNNs of \(2\times 20\) hidden layers} \\ \hline \(\aleph\) & OSIF Model & \(k\) & Acc. & OSIF Model & \(k\) & Acc. \\ \hline \multirow{2}{*}{50} & \(2\times 10\) & 1 & 76.6\% & \(2\times 20\) & 1 & 77.2\% \\  & \(2\times 10\) & 2 & **77.0\%** & \(2\times 20\) & 2 & **77.9\%** \\ \hline \multirow{2}{*}{100} & \(2\times 10\) & 1 & 81.1\% & \(2\times 20\) & 1 & 81.3\% \\  & \(2\times 10\) & 2 & **82.3\%** & \(2\times 20\) & 2 & 81.8\% \\  & \(2\times 10\) & 3 & 82.2\% & \(2\times 20\) & 3 & **82.1\%** \\ \hline \multicolumn{8}{c}{DNNs of \(2\times 40\) hidden layers} \\ \hline \(\aleph\) & OSIF Model & \(k\) & Acc. & OSIF Model & \(k\) & Acc. \\ \hline \multirow{2}{*}{50} & \(2\times 10\) & 1 & 78.3\% & \(2\times 20\) & 1 & 78.6\% \\  & \(2\times 10\) & 2 & **79.4\%** & \(2\times 20\) & 2 & **79.6\%** \\ \hline \multirow{2}{*}{100} & \(2\times 10\) & 1 & 82.4\% & \(2\times 20\) & 1 & 82.6\% \\  & \(2\times 10\) & 2 & 83.1\% & \(2\times 20\) & 2 & 83.2\% \\ \multicolumn{8}{c}{\(2\times 10\)} & \(3\) & **83.7\%** & \(2\times 20\) & 3 & **84.0\%** \\ \hline \end{tabular}
\end{table}
Table 2: Accuracy of DNNs of different architectures with sparse input features selected by Alg. 1 on FashionMNIST
#### 3.1.2 Comparisons between feature selection approaches

In this section, we compare the performance of the MILP-based features-selection approach (_i.e.,_ Alg. 1) to some other simple feature-selection approaches. The other feature selection techniques considered in the comparison are random feature selection, data-based feature selection, and DNN weights-based. In the following paragraphs, we briefly describe the feature selection algorithms that are used as reference points for the comparison.

The random feature selection uniformly samples a subset of \(\aleph\) input features, and we present the average accuracy with the standard deviation over 5 DNNs trained with input features independently selected by this approach. This approach is included, as it is the simplest approach to down-sample the full input.

The data-based feature selection is conducted in the following way: i) we first calculate the mean value of each input feature over the whole train dataset; ii) the features with the largest \(\aleph\) mean are selected. The motivation behind this simple heuristic is that it selects the pixels that are most strongly colored over the training data, _i.e._, the strongest signals. For example, selecting a pixel uncolored in all images of the training data does not make sense as that input does not contain any information for the training data.

In the DNN weights-based approach, we use the same DNN models as we use in the MILP-based selection, but the inputs are now selected based on the weights of the inputs. For each input, we sum up the absolute values of all the weights from the input to the nodes in the consecutive layer and select the ones with the largest sum. This can be seen as a form of pruning of inputs, and the motivation is that inputs with small, or almost zero, weights should be less important as these inputs have less impact in the DNN.

In Table 3, we compare MILP-based feature selection (_i.e.,_ Alg. 1) with random selection, data-based selection, and weights-based selection. The figure for Alg. 1 has the best result from Table 1, where the accuracy of DNNs with sparse input features is only 5 points with \(\aleph=100\) less than the accuracy with full input features. It can be observed that our method has the best overall performance. For MNIST, the random selection has the worst

Figure 3: MNIST and FashionMNIST input features selected by Alg. 1 with \(\aleph=100\) and \(k\in\{1,2,3\}\)
performance, but the data-based selection and weights-based selection achieve a slightly worse performance than our method.

Table 4 compares the performance of features selections on FashionMNIST. The results show a different pattern to MNIST. Our methods still have the best overall performance over different settings by maintaining the accuracy to 84% (resp. 88%) with \(\aleph=100\) for DNNs \(2\times 20\) (resp. \(2\times 40\)), while the accuracy of DNNs with full input features are 86% (resp. 88%). While delivering the worst performance on MNIST, random selection has a very close performance to our method on FashionMNIST. The weights-based selection still lies in third overall, while the data-based selection is much worse than the other three methods (_e.g.,_ 58% for the DNN \(2\times 20\) with \(\aleph=100\) and 62% for the DNN \(2\times 40\) with \(\aleph=100\)).

The weights-based selection performs decently on both data sets compared to random and data-based selection. However, based on the results it is clear that MILP-based selection (_i.e.,_ Alg. 1) can extract more knowledge from the DNN model regarding the importance of inputs compared to simply analyzing the weights. The overall performance of MILP-based selection is more stable than other feature selection methods on both datasets.

### Robustness to adversarial inputs

The robustness of a trained DNN classifier can also be analyzed using MILP, _e.g.,_ in verification or finding adversarial input. We use the minimal distorted adversary as a measure of model robustness \(x\) under \(l_{\infty}\) norm [10]. For a given image \(x_{\text{image}}\), the minimal adversary problem [29] can be formulated as

\[\begin{split}\min\ \epsilon&\\ \text{s.t.(\ref{eq:x_1})-(\ref{eq:x_2})},x_{i}^{I}\leq x_{j}^{L}, ||x^{1}-x_{\text{image}}||_{\infty}\leq\epsilon,\end{split} \tag{7}\]

\begin{table}
\begin{tabular}{c|c c c c|c} \hline \multirow{2}{*}{DNN} & \multicolumn{5}{c|}{Feature Selection Approaches} & \multirow{2}{*}{\(\aleph=784\)} \\ \cline{2-2} \cline{4-6}  & \(\aleph\) & MILP-based & & Random & Data-based & Weights-based \\ \hline \multirow{2}{*}{\(2\times 20\)} & 50 & 77.9\% & \(77.2\pm 0.8\%\) & 49.6\% & 73.8\% & \multirow{2}{*}{86.3\%} \\  & 100 & 82.3\% & \(80.6\pm 1.1\%\) & 58.4\% & 80.3\% & \\ \hline \multirow{2}{*}{\(2\times 40\)} & 50 & 79.6\% & \(78.8\pm 0.4\%\) & 51.6\% & 74.6\% & \multirow{2}{*}{87.5\%} \\  & 100 & 84.0\% & \(82.8\pm 0.4\%\) & 62.3\% & 81.7\% & \\ \hline \end{tabular}
\end{table}
Table 4: Accuracy of DNNs with sparse input features selected by different methods on FashionMNIST

\begin{table}
\begin{tabular}{c|c c c c c|c} \hline \multirow{2}{*}{DNN} & \multicolumn{5}{c|}{Feature Selection Approaches} & \multirow{2}{*}{\(\aleph=784\)} \\ \cline{2-2} \cline{4-6}  & \(\aleph\) & MILP-based & & Random & Data-based & Weights-based \\ \hline \multirow{2}{*}{\(2\times 20\)} & 50 & 86.6\% & \(77.4\pm 2.2\%\) & 81.3\% & 80.7\% & \multirow{2}{*}{95.7\%} \\  & 100 & 91.2\% & \(86.0\pm 2.5\%\) & 89.2\% & 89.0\% & \\ \hline \multirow{2}{*}{\(2\times 40\)} & 50 & 89.2\% & \(76.0\pm 3.9\%\) & 85.1\% & 83.3\% & \multirow{2}{*}{97.1\%} \\  & 100 & 93.4\% & \(90.6\pm 1.1\%\) & 92.4\% & 91.2\% & \\ \hline \end{tabular}
\end{table}
Table 3: Accuracy of DNNs with sparse input features selected by different methods on MNIST
where \(i\) is the true label of image \(x_{\text{image}}\) and \(j\) is an adversarial label. Simply put, problem (7) finds the smallest perturbation, defined by the \(\ell_{\infty}\) norm, such that the trained DNN erroneously classifies image \(x_{\text{image}}\) as the adversarial label \(j\). We hypothesize that DNNs trained with fewer (well-selected) features are more robust to such attacks, as there are fewer inputs as degrees of freedom. Furthermore, we note that the robustness of smaller DNNs can be analyzed with significantly less computational effort.

Table 5 shows the minimal adversarial distance (mean and standard deviation over 100 instances), defined by (7) for DNNs trained on MNIST and FashionMNIST with MILP-based feature selection. The adversaries are generated for the first 100 instances of the respective test datasets, with adversarial labels selected randomly. Furthermore, we report the mean percentage increase, \(\Delta\), in minimal adversarial distance over the 100 instances for the reduced input DNNs compared to the full input. In all cases, reducing the the number of inputs \(\aleph\) results in a more robust classifier. For the \(2\times 40\) DNN trained on FashionMNIST, reducing the number of inputs from 784 to 50 increases the mean minimal adversarial distance by almost 90%, with a loss in accuracy of \(<\)10%.

## 4 Conclusion

In the paper, we have presented an MILP-based framework using trained DNNs to extract information about salient input features. The proposed algorithm is able to drastically reduce the size of the input by using the input features that are most important for each category according to the DNN, given a regularization on the input size and spread of selected features. The numerical results show that the proposed algorithm is able to efficiently select a small set of features for which a good prediction accuracy can be obtained. The results also show that the proposed input feature selection can improve the robustness toward adversarial attacks.
