# A Robust Optimization Approach to Deep Learning

Dimitris Bertsimas<EMAIL_ADDRESS>
Sloan School of Management and Operations Research Center
Massachusetts Institute of Technology
Cambridge, MA 02139, USA Xavier Boix<EMAIL_ADDRESS>
Department of Brain and Cognitive Sciences
Massachusetts Institute of Technology
Cambridge, MA 02139, USA Kimberly Villalobos Carballo<EMAIL_ADDRESS>
Operations Research Center
Massachusetts Institute of Technology
Cambridge, MA 02139, USA Dick den Hertog<EMAIL_ADDRESS>
Amsterdam Business School
University of Amsterdam
Plantage Muidergracht 12, 1018 TV Amsterdam, The Netherlands

###### Abstract

Many state-of-the-art adversarial training methods leverage upper bounds of
the adversarial loss to provide security guarantees. Yet, these methods
require computations at each training step that can not be incorporated in the
gradient for backpropagation. We introduce a new, more principled approach to
adversarial training based on a closed form solution of an upper bound of the
adversarial loss, which can be effectively trained with backpropagation. This
bound is facilitated by state-of-the-art tools from robust optimization. We
derive two new methods with our approach. The first method (Approximated
Robust Upper Bound or aRUB) uses the first order approximation of the network
as well as basic tools from linear robust optimization to obtain an
approximate upper bound of the adversarial loss that can be easily
implemented. The second method (Robust Upper Bound or RUB), computes an exact
upper bound of the adversarial loss. Across a variety of tabular and vision
data sets we demonstrate the effectiveness of our more principled approach
—RUB is substantially more robust than state-of-the-art methods for larger
perturbations, while aRUB matches the performance of state-of-the-art methods
for small perturbations. Also, both RUB and aRUB run faster than standard
adversarial training (at the expense of an increase in memory). All the code
to reproduce the results can be found at https://github.com/kimvc7/Robustness.

Keywords: deep learning, optimization, robustness

##  1 Introduction

Robustness of neural networks for classification problems has received
increasing attention in the past few years, since it was exposed that these
models could be easily fooled by introducing some small perturbation in the
input data. These perturbed inputs, which are commonly referred to as
adversarial examples, are visually indistinguishable from the natural input,
and neural networks simply trained to maximize accuracy often assign them to
an incorrect class (Szegedy et al., 2013). This problem becomes particularly
relevant when considering applications related to self-driving cars or
medicine, in which adversarial examples represent an important security threat
(Kurakin et al., 2016).

The Machine Learning community has recently developed multiple heuristics to
make neural networks robust. The most popular ones are perhaps those based on
adversarial training, a method first proposed by Goodfellow et al. (2015) and
which consists in training the neural network using adversarial inputs instead
of the standard data. The defense by Madry et al. (2019), which finds the
adversarial examples with bounded norm using iterative projected gradient
descent (PGD) with random starting points, has proved to be one of the most
effective methods (Tjeng et al., 2019). Other heuristic defenses rely on
preprocessing or projecting the input space (Lamb et al., 2018; Kabilan et
al., 2018; Ilyas et al., 2017), on randomizing the neurons (Prakash et al.,
2018; Xie et al., 2017) or on adding a regularization term to the objective
function (Ross and Doshi-Velez, 2017; Hein and Andriushchenko, 2017; Yan et
al., 2018). There is a plethora of heuristics for adversarial robustness by
now. Yet, these defenses are only effective to adversarial attacks of small
magnitude and are vulnerable to attacks of larger magnitude or to new attacks
(Athalye et al., 2018).

A recent strand of research has been to leverage upper bounds of the
adversarial loss to improve adversarial robustness. These upper bounds provide
security guarantees against adversarial attacks, even new ones, by finding a
mathematical proof that a network is not susceptible to any attack, e.g.
(Dathathri et al., 2020; Raghunathan et al., 2018b; Katz et al., 2017; Tjeng
et al., 2019; Bunel et al., 2017; Anderson et al., 2020; Singh et al., 2018;
Zhang et al., 2018; Weng et al., 2018; Gehr et al., 2018; Dvijotham et al.,
2018b; Lecuyer et al., 2019; Cohen et al., 2019). The upper bounds have also
been incorporated in the training of adversarial defenses. Wong and Kolter
(2018) for instance, find an upper bound for the adversarial loss by
considering the polytope generated by adversarial examples with bounded norm
and minimizing the loss over a convex polytope that contains it. An upper
bound on the adversarial loss is also computed in Raghunathan et al. (2018a)
by solving instead a semidefinite program. Other more scalable and effective
methods based on minimizing an upper bound of the adversarial loss have been
introduced (Gowal et al., 2019; Balunovic and Vechev, 2019; Mirman et al.,
2018; Dvijotham et al., 2018a; Wong et al., 2018; Zhang et al., 2019).

While these methods provide security guarantees, most of them require complex
intermediate computations that depend on the training parameters of the
network and can not be incorporated in the gradient for backpropagation. Since
the direction of the gradient is then not precise, the optimization process
may be problematic. Even though a couple of these approaches do provide an
upper bound with closed form solution (Raghunathan et al., 2018a; Gowal et
al., 2019), they only work for specific cases (Raghunathan et al. (2018a) only
works for neural networks with two layers, and the upper bound in Gowal et al.
(2019) is often loose and it only applies when the perturbations are in the
$L_{\infty}$ norm-bounded ball). Given this state of affairs, we believe we
need a new, principled approach to optimize an upper bound of the adversarial
loss.

A promising approach that has not been comprehensively explored is the usage
of state-of-the-art robust optimization tools (Bertsimas and den Hertog,
2022). Robust optimization has been tremendously effective to solve
optimization problems with uncertainty in the parameters, which may arise from
rounding or implementation errors. It provides guarantees to uncertainty in
the parameters of the model and could be leveraged for deep learning
(Bertsimas et al., 2019). In this paper, we use state-of-the-art robust
optimization tools to derive a new closed form solution of an upper bound of
the adversarial loss. This leads to a more principled approach to adversarial
training that can be effectively trained with backpropagation.

We develop two new methods for training deep learning models that are robust
against input perturbations. The first method (Approximated Robust Upper Bound
or aRUB), minimizes an approximated upper bound of the robust loss for $L_{p}$
norm bounded uncertainty sets, for general $p$. It is simple to implement and
performs similarly to adversarial training on small uncertainty sets while
running significantly faster. The second method (Robust Upper Bound or RUB),
minimizes an exact upper bound of the robust loss specifically for $L_{1}$
norm bounded uncertainty sets. This method shows the best performance for
larger uncertainty sets and also runs faster than adversarial training. More
importantly, it provides security guarantees against $L_{1}$ norm bounded
adversarial attacks. Lastly, we show that robust training can significantly
improve the natural accuracy compared to standard non-robust training methods.
More concretely, we introduce the following robustness methods:

  * •

Approximated Robust Upper Bound or aRUB: We develop a simple method to
approximate an upper bound of the robust problem by adding a regularization
term for each class separately. As an alternative to standard adversarial
training (which relies on linear approximations to find good adversarial
attacks), we use the first order approximation of the network to estimate the
worst case scenario for each individual class. We then apply standard results
from linear robust optimization to obtain a new objective that behaves like an
upper bound of the robust loss and which can be tractably minimized for robust
training. This method can be easily implemented and performs very well when
the uncertainty set radius $\rho$ is small. It has a higher memory requirement
than nominal training, but it significantly reduces training time compared to
other robust methods like adversarial training.

  * •

Robust Upper Bound or RUB: We extended state-of-the-art tools from robust
optimization to function that like neural networks are neither convex nor
concave. By splitting each layer of the network as the sum of a convex
function and a concave function, we are able to obtain an exact upper bound of
the robust loss for the case in which the uncertainty set is the $L_{1}$
sphere. Since the dual function of the $L_{1}$ norm is the $L_{\infty}$ norm,
we convert the maximum over the uncertainty set into a maximum over a finite
set. In the end, instead of minimizing the worst case loss over an infinite
uncertainty set, the new objective minimizes the worst case loss over a
discrete set whose cardinality is twice the dimension of the input data. While
this represents a significant increase in memory for high dimensional inputs,
we show that this approach remains tractable for multiple applications. In
addition, this methods allows for parallelization and the overall
computational time is much lower than the one achieved with other algorithms
like adversarial training. The main advantage of this method is that it
provides security guarantees against adversarial examples bounded in the
$L_{1}$ norm. Additionally, we also show experimentally that this method
generally achieves the highest adversarial accuracies for larger uncertainty
sets.

Also, we show that these methods consistently achieve higher standard accuracy
(i.e., non adversarial accuracy), than the nominal neural networks trained
without robustness. While this result is not true for a general choice of
uncertainty set (see for example Ilyas et al. (2019)), we observe that when
the uncertainty set has the appropriate size it can significantly improve the
classification performance of the network, which is consistent with the
results obtained for other classification models like Support Vector Machines,
Logistic Regression and Classification Trees (Bertsimas et al., 2019).

The paper is organized as follows: Section 2 revisits related work with robust
Optimization and deep learning, Section 3 defines the robust problem, Section
4 presents our first method (Approximate Robust Upper Bound), and Section 5
contains our second method (Robust Upper Bound). Lastly, Section 6 contains
the results for the computational experiments.

##  2 Previous Works on Robust Optimization

Over the last two decades, robust optimization has become a successful
approach to solve optimization problems under uncertainty. For an overview of
the primary research in this field we refer the reader to Bertsimas et al.
(2011). Areas like mathematical programming and engineering have long applied
these tools to develop models that are robust against uncertainty in the
parameters, which may arise from rounding or implementation errors. For many
applications, the robust problem can be reformulated as a tractable
optimization problem, which is referred to as the robust counterpart. For
instance, for several types of uncertainty sets, the robust counterpart of a
linear programming problem can be written as a linear or conic programming
problem (Ben-Tal et al., 2009), which can be solved with many of the current
optimization software. While there is not a systematic way to find robust
counterparts for a general nonlinear uncertain problem, multiple techniques
have been developed to obtain tractable formulations in some specific
nonlinear cases.

As shown in Ben-Tal et al. (2009), the exact robust counterpart is known for
Conic Quadratic Problems and Semidefinite problems in which the uncertainty is
finite, an interval or an unstructured norm-bounded set. More generally, it is
shown in Ben-Tal et al. (2015) that for problems in which the objective
function or the constraints are concave in the uncertainty parameters, Fenchel
duality can be used to exactly derive the corresponding robust counterpart.
While the result does not necessarily have a closed form, the authors show
that it yields a tractable formulation for the most common uncertainty sets
(e.g. polyhedral and ellipsoidal uncertainty sets).

The problem becomes significantly more complex when the functions in the
objective or in the constraints are instead convex in the uncertainty
(Chassein and Goerigk, 2019). Since obtaining exact robust counterparts in
these cases is generally infeasible, safe approximations are considered
instead (Bertsimas et al., 2020). For instance, Zhen et al. (2017) develop
safe approximations for the specific cases of second order cone and
semidefinite programming constraints with polyhedral uncertainty. These
techniques are generalized in Roos et al. (2020), where the authors convert
the robust counterpart to an adjustable robust optimization problem that
produces a safe approximation for any problem that is convex in the
optimization variables as well as in the the uncertain parameters.

Even though the approaches mentioned above consider uncertainty in the
parameters of the model as opposed to uncertainty in the input data, the same
techniques can be utilized for obtaining robust counterparts in the latter
case. In fact, the robust optimization methodologies have recently been
applied to develop machine learning models that are robust against
perturbations in the input data. In Bertsimas et al. (2019), for example, the
authors consider uncertainty in the data features as well as in the data
labels to obtain tractable robust counterparts for some of the major
classification methods: support vector machines, logistic regression, and
decision trees. However, due to the high complexity of neural networks as well
as the large dimensions of the problems in which they are often utilized,
robust counterparts or safe approximations for this type of models have not
yet been developed. There are two major challenges with applying robust
optimization tools for training robust neural networks:

  1. (i)

Neural networks are neither convex nor concave functions: As mentioned
earlier, robust counterparts are difficult to find for a general problem.
Although plenty of work has been done to find tractable reformulations as well
as safe approximations, all of them rely on the underlying function being a
convex or a concave function of the uncertainty parameters. Unfortunately,
neural networks don’t satisfy either condition, which makes it really
difficult to apply any of the approaches discussed above.

  2. (ii)

The robust counterpart needs to preserve the scalability of the model: Neural
networks are most successful in problems involving vision data sets, which
often imply large input dimensions and enormous amount of data. For the most
part, they can still be successfully trained thanks to the fact that back
propagation algorithms can be applied to solve the corresponding unconstrained
optimization problem. However, the robust optimization techniques for both
convex and concave cases often require the addition of new constraints and
variables for each data sample, increasing significantly the number of
parameters of the network and making it very difficult to use standard machine
learning software for training.

A straightforward way to overcome both of these difficulties would be to
replace the loss function of the network with its first order approximation.
However, this loss function is usually highly nonlinear and therefore the
linear approximation is very inaccurate. Our method aRUB explores a slight
modification of this approach that significantly improves adversarial accuracy
by considering only the linear approximation of the network’s output layer.

Alternatively, a more rigorous approach to overcome problem (i) would be to
piece-wise analyze the convexity of the network and apply the robust
optimization techniques in each piece separately, but this approach would
introduce additional variables that are in conflict with requirement (ii). For
our RUB method we then develop a general framework to split the network by
convexity type, and we show that in the specific case in which the uncertainty
set is the $L_{1}$ norm bounded sphere, we can solve for the extra variables
and obtain an unconstrained problem that can be tractably solved using
standard gradient descent techniques.

##  3 The Robust Optimization Problem

We consider a classification problem over data points
${\boldsymbol{x}}\in\mathbb{R}^{M}$ with $K$ different classes. Given weight
matrices ${\boldsymbol{W}}^{\ell}\in\mathbb{R}^{r_{\ell-1}\times r_{\ell}}$
and bias vectors ${\boldsymbol{b}}^{\ell}\in\mathbb{R}^{r_{\ell}}$ for
$\ell\in[L]$, such that $r_{0}=M,r_{L}=K$, the corresponding feed forward
neural network with $L$ layers and ReLU activation function is defined by the
equations

$\displaystyle{\boldsymbol{z}}^{1}(\theta,{\boldsymbol{x}})$
$\displaystyle={\boldsymbol{W}}^{1}{\boldsymbol{x}}+{\boldsymbol{b}}^{1},$
$\displaystyle{\boldsymbol{z}}^{\ell}(\theta,{\boldsymbol{x}})$
$\displaystyle={\boldsymbol{W}}^{\ell}[{\boldsymbol{z}}^{\ell-1}(\theta,{\boldsymbol{x}})]^{+}+{\boldsymbol{b}}^{\ell},\quad\forall\>2\leq\ell\leq
L,$

where $\theta$ denotes the set of parameters
$({\boldsymbol{W}}^{\ell},{\boldsymbol{b}}^{\ell})$ for all $\ell\in[L]$ and
$[{\boldsymbol{x}}]^{+}$ is the result of applying the ReLU function to each
coordinate of ${\boldsymbol{x}}$. For fixed parameters $\theta$, the network
assigns a sample ${\boldsymbol{x}}$ to the class
$\hat{y}=\arg\max_{k}\>{\boldsymbol{z}}^{L}_{k}(\theta,{\boldsymbol{x}})$. And
given a data set $\\{({{\boldsymbol{x}}}_{n},{{y}}_{n})\\}_{n=1}^{N}$, where
$y_{n}$ is the target class of ${\boldsymbol{x}}_{n}$, the optimal parameters
$\theta$ are usually found by minimizing the empirical loss

$\displaystyle\min_{\theta}\enspace\frac{1}{N}{\sum_{n=1}^{N}\mathcal{L}({y_{n}},{\boldsymbol{z}}^{L}(\theta,{{\boldsymbol{x}}}_{n}))},$
(1)

with respect to a specific loss function $\mathcal{L}$.

In the robust optimization framework, however, we want to find the parameters
$\theta$ by minimizing the worst case loss achieved over an uncertainty set of
the input. More specifically, instead of solving the nominal problem in Eq.
(1), we want to solve the following min-max problem:

$\displaystyle\min_{\theta}\enspace\frac{1}{N}{\sum_{n=1}^{N}\max_{\boldsymbol{\delta}\in\mathcal{U}}\>\mathcal{L}\left({y_{n}},{\boldsymbol{z}}^{L}(\theta,\>{{\boldsymbol{x}}}_{n}+\boldsymbol{\delta})\right)},$
(2)

for some uncertainty set $\mathcal{U}\subset\mathbb{R}^{M}$.
$\in\mathbb{R}^{K}$

Unfortunately, a closed form expression for the inner maximization problem
above is unknown and solving the min-max problem is notoriously difficult.
Robust optimization provides multiple tools for solving such problems when the
loss function is either convex or concave in the input variables. For example,
in the case of concave loss functions, a common approach would be to take the
dual of the maximization problem so that the problem can be formulated as a
single minimization problem (Bertsimas and den Hertog, 2022). If the loss
function is instead convex, Fenchel’s duality as well as conjugate functions
can be used to find upper bounds and lower bounds of the maximization problem.
However, there is no general framework developed for loss functions that do
not fall into those categories, like in the case of neural networks. In this
paper, we will focus on the specific case in which the uncertainty set is the
ball of radius $\rho$ in the $\text{L}_{p}$ space; i.e.
$\mathcal{U}=\\{{\boldsymbol{\delta}}:\|{\boldsymbol{\delta}}\|_{p}\leq\rho\\}$.

##  4 Approximate Robust Upper Bound for small $\rho$

Perhaps the most intuitive approach to tackle problem (2) is to consider the
first order approximation of the loss function

$\displaystyle{\mathcal{L}}\left(y,{\boldsymbol{z}}^{L}(\theta,{{\boldsymbol{x}}}+\boldsymbol{\delta})\right)\approx\mathcal{L}\left(y,{\boldsymbol{z}}^{L}(\theta,\>{{\boldsymbol{x}}})\right)+\boldsymbol{\delta}^{\top}\nabla_{{\boldsymbol{x}}}\mathcal{L}\left(y,{\boldsymbol{z}}^{L}(\theta,\>{{\boldsymbol{x}}})\right),$

since the right hand side is a linear function of $\boldsymbol{\delta}$ and
the maximization problem can be more easily solved for linear functions of the
uncertainty. For example, it is not hard to see that the first order
approximation reaches its maximum value in
$\mathcal{U}=\\{\boldsymbol{\delta}:\|\boldsymbol{\delta}\|_{\infty}\leq\rho\\}$
exactly at
$\boldsymbol{\delta}^{\star}=\rho\>sign\left(\nabla_{{\boldsymbol{x}}}\mathcal{L}\big{(}y,{\boldsymbol{z}}^{L}(\theta\>{{\boldsymbol{x}}})\big{)}\right)$.
This approach is referred to as fast gradient sign method and it was first
explored in Goodfellow et al. (2015), where the networks are trained with
adversarial examples generated as
$\boldsymbol{x}+\boldsymbol{\delta}^{\star}$. A similar approach was proposed
in Huang et al. (2016), where the authors considered the cross entropy loss
and use the linear approximation of the softmax layer instead of the
approximation of the entire loss function. In these examples, linear
approximations are used to find near optimal perturbations that can produce
strong adversarial examples for training, but not to modify the objective
function.

An alternative to these methods would then be to train the network with the
natural data and replace the loss function with its linear approximation,
transforming the problem into

$\displaystyle\min_{\theta}\enspace\frac{1}{N}{\sum_{n=1}^{N}\>\max_{\boldsymbol{\delta}\in\mathcal{U}}\>\mathcal{L}\left(y_{n},{\boldsymbol{z}}^{L}(\theta,{{\boldsymbol{x}}}_{n})\right)+\boldsymbol{\delta}^{\top}\nabla_{{\boldsymbol{x}}}\mathcal{L}\left(y_{n},{\boldsymbol{z}}^{L}(\theta,{{\boldsymbol{x}}}_{n})\right)}$
$\displaystyle=$
$\displaystyle\min_{\theta}\enspace\frac{1}{N}{\sum_{n=1}^{N}\>\mathcal{L}\left(y_{n},{\boldsymbol{z}}^{L}(\theta,{{\boldsymbol{x}}}_{n})\right)+\rho\|\nabla_{{\boldsymbol{x}}}\mathcal{L}\left(y_{n},{\boldsymbol{z}}^{L}(\theta,{{\boldsymbol{x}}}_{n})\right)\|_{q}},$
(3)

where $\|\>\|_{q}$ is the dual norm of $\|\>\|_{p}$, satisfying
$\frac{1}{p}+\frac{1}{q}=1$. However, this approach has not been very popular,
since the loss function is highly nonlinear and the adversarial accuracy
obtained by optimizing the linear approximation is much worse than the one
obtained with adversarial training even for small values of $\rho$ (see the
Baseline method in Section 6). Nonetheless, from the network definitions in
Section 3 we can see that each component of the network
${\boldsymbol{z}}^{L}(\theta,\boldsymbol{x})$ is in fact a continuous
piecewise linear function, suggesting that the first order approximation of
${\boldsymbol{z}}^{L}$ is more precise than that of
$\mathcal{L}(y,{\boldsymbol{z}}^{L})$. In fact, we expect the outputs
${\boldsymbol{z}}^{L}(\theta,\boldsymbol{x})$ and
${\boldsymbol{z}}^{L}(\theta,\boldsymbol{x}+\boldsymbol{\delta})$ to be in the
same linear piece when $\boldsymbol{x}+\boldsymbol{\delta}$ is close to
$\boldsymbol{x}$. In other words, the linear approximation

$\displaystyle{\boldsymbol{z}}^{L}(\theta,\boldsymbol{x}+\boldsymbol{\delta})\approx{\boldsymbol{z}}^{L}(\theta,\boldsymbol{x})+\boldsymbol{\delta}^{\top}\nabla_{{\boldsymbol{x}}}{\boldsymbol{z}}^{L}(\theta,\boldsymbol{x})$
(4)

is exact for small enough $\boldsymbol{\delta}$. We can then approximately
solve the adversarial problem for each class $k$ as

$\displaystyle\max_{\boldsymbol{\delta}\in\mathcal{U}}\>({\boldsymbol{e}}_{k}-{\boldsymbol{e}}_{y})^{\top}{\boldsymbol{z}}^{L}(\theta,\boldsymbol{x}+\boldsymbol{\delta})$
$\displaystyle\approx\max_{\boldsymbol{\delta}\in\mathcal{U}}\>({\boldsymbol{e}}_{k}-{\boldsymbol{e}}_{y})^{\top}{\boldsymbol{z}}^{L}(\theta,\boldsymbol{x})+\boldsymbol{\delta}^{\top}\nabla_{x}({\boldsymbol{e}}_{k}-{\boldsymbol{e}}_{y})^{\top}{\boldsymbol{z}}^{L}(\theta,\boldsymbol{x})$
$\displaystyle=({\boldsymbol{e}}_{k}-{\boldsymbol{e}}_{y})^{\top}{\boldsymbol{z}}^{L}(\theta,\boldsymbol{x})+\rho\|({\boldsymbol{e}}_{k}-{\boldsymbol{e}}_{y})^{\top}\nabla_{x}{\boldsymbol{z}}^{L}(\theta,\boldsymbol{x})\|_{q}$
$\displaystyle=\boldsymbol{c}_{k}^{\top}{\boldsymbol{z}}^{L}(\theta,\boldsymbol{x})+\rho\|{\boldsymbol{c}}_{k}^{\top}\nabla_{x}{\boldsymbol{z}}^{L}(\theta,\boldsymbol{x})\|_{q},$
(5)

where $\boldsymbol{e}_{k}$ refers to the one-hot vectors with a $1$ in the
$k^{th}$ coordinate and $0$ everywhere else, and $\boldsymbol{c}_{k}$ is
defined as $\boldsymbol{c}_{k}\coloneqq\boldsymbol{e}_{k}-\boldsymbol{e}_{y}$.
Since the loss functions utilized in deep learning for classification problems
(cross-entropy loss, hinge loss, and zero-one loss) are monotone and
transitionally invariant; i.e. they satisfy

$\displaystyle\mathcal{L}(y,{\boldsymbol{z}})=\mathcal{L}(y,{\boldsymbol{z}}-c\boldsymbol{1})\quad\forall\>c\in\mathbb{R},$
(6)

we can use Eq. (5) to obtain an approximate upper bound of problem (2) as
follows:

$\displaystyle\min_{\theta}\>\max_{\boldsymbol{\delta}\in\mathcal{U}}\>\mathcal{L}(y,{\boldsymbol{z}}^{L}(\theta,{\boldsymbol{x}}+\boldsymbol{\delta}))$
(7) $\displaystyle=$
$\displaystyle\min_{\theta}\>\max_{\boldsymbol{\delta}\in\mathcal{U}}\>\mathcal{L}\left(y,{\boldsymbol{z}}^{L}(\theta,{\boldsymbol{x}}+\boldsymbol{\delta})-{\boldsymbol{e}}_{y}^{\top}{\boldsymbol{z}}^{L}(\theta,{\boldsymbol{x}}+\boldsymbol{\delta})\boldsymbol{1}\right),\hskip
128.0374pt(\text{By Eq. \eqref{trans-invariant}})$ (8) $\displaystyle\leq$
$\displaystyle\min_{\theta}\>\mathcal{L}\left(y,\bigg{(}\max_{\boldsymbol{\delta}\in\mathcal{U}}{\boldsymbol{c}}_{1}^{\top}{\boldsymbol{z}}^{L}(\theta,{\boldsymbol{x}}+\boldsymbol{\delta}),\dots,\max_{\boldsymbol{\delta}\in\mathcal{U}}{\boldsymbol{c}}_{K}^{\top}{\boldsymbol{z}}^{L}(\theta,{\boldsymbol{x}}+\boldsymbol{\delta}\bigg{)}\right),\quad\hskip
54.06006pt(\text{By monoticity})$ (9) $\displaystyle\approx$
$\displaystyle\min_{\theta}\>\mathcal{L}\left(y,\big{(}\boldsymbol{c}_{1}^{\top}{\boldsymbol{z}}^{L}(\theta,{\boldsymbol{x}})+\rho\|\boldsymbol{c}_{1}^{\top}\nabla_{x}{\boldsymbol{z}}^{L}(\theta,{\boldsymbol{x}})\|_{q},\dots,\boldsymbol{c}_{K}^{\top}{\boldsymbol{z}}^{L}(\theta,{\boldsymbol{x}})+\rho\|\boldsymbol{c}_{K}^{\top}\nabla_{x}{\boldsymbol{z}}^{L}(\theta,{\boldsymbol{x}})\|_{q}\big{)}\right).$
(10) $\displaystyle\hskip 361.3499pt(\text{By Eq. \eqref{class-adversary}})$

Therefore, we propose to train the network by minimizing Eq. (10) instead of
the standard average loss. For the particular case of the cross entropy loss
with softmax activation function in the output layer, the exact optimization
problem to be solved would be

$\displaystyle\min_{\theta}\>\frac{1}{N}\sum_{n=1}^{N}\log\left(\sum_{k}e^{(\boldsymbol{e}_{k}-\boldsymbol{e}_{y_{n}})^{\top}{\boldsymbol{z}}^{L}(\theta,{\boldsymbol{x}})+\rho\|\nabla_{{\boldsymbol{x}}}(\boldsymbol{e}_{k}-\boldsymbol{e}_{y_{n}})^{\top}{\boldsymbol{z}}^{L}(\theta,{\boldsymbol{x}})\|_{q}}\right).$
(11)

Even though this expression is not guaranteed to be an upper bound of Eq. (2),
we observe across a variety of experiments that Eq. (11) behaves like an upper
bound of the adversarial loss (see Figure 1 and Table 1), indicating that the
upper bound obtained by monotonocity in Eq. (9) generally counteracts the
error from the first order approximation of ${\boldsymbol{z}}^{L}$.

(a) MNIST, $\rho=0.01$

(b) CIFAR, $\rho=0.001$

(c) Fashion MNIST, $\rho=0.01$

(d) MNIST, $\rho=0.1$

(e) CIFAR, $\rho=0.01$

(f) Fashion MNIST, $\rho=0.1$

(g) MNIST, $\rho=0.3$

(h) CIFAR, $\rho=0.03$

(i) Fashion MNIST, $\rho=0.3$

Figure 1: Adversarial loss (cross entropy loss evaluated at adversarial images bounded in $L_{\infty}$ norm) vs the approximate upper bound (Eq. (11) evaluated at the natural images). The columns correspond to experiments with the MNIST, CIFAR, and Fashion MNIST data sets from left to right, and the value of $\rho$ increases for lower rows. |  $\boldsymbol{\rho}$ = $\boldsymbol{0.0}$ | $\boldsymbol{0.0008}$ | $\boldsymbol{0.001}$ | $\boldsymbol{0.0015}$ | $\boldsymbol{0.002}$ | $\boldsymbol{0.003}$ | $\boldsymbol{0.01}$ | $\boldsymbol{0.1}$ | $\boldsymbol{0.3}$ | $\boldsymbol{0.5}$ | $\boldsymbol{1.0}$
---|---|---|---|---|---|---|---|---|---|---|---
$L_{\infty}$ | 0.94 | 0.99 | 0.99 | 0.99 | 0.99 | 0.99 | 0.99 | 0.95 | 0.86 | 0.81 | 0.79
$L_{1}$ | 0.933 | 0.99 | 0.99 | 0.99 | 0.99 | 0.99 | 0.99 | 0.99 | 0.99 | 0.99 | 0.98
Table 1: Percentage of times when the aRUB approach yields an upper bound of
the adversarial cross entropy with respect to PGD attacks. For each row, both
training methods use the $L_{p}$ norm indicated on the first column.
Percentages are computed across multiple networks trained with 46 UCI data
sets.

##  5 Robust Upper bound for the $L_{1}$ norm and general $\rho$.

In this section we derive an exact upper bound for the robust counterpart of
the inner maximization problem in Eq. (2) for the specific case in which the
uncertainty set is the ball of radius $\rho$ in the $\text{L}_{1}$ space; i.e.
$\mathcal{U}=\\{{\boldsymbol{\delta}}:\|{\boldsymbol{\delta}}\|_{1}\leq\rho\\}$.
To find these upper bound we will use several tools from nonlinear and robust
optimization. For this, we will make use of three lemmas that rely on the
definition of conjugate functions as well as the Fenchel duality theorem, and
their proofs can be found in Appendix A. Together, these lemmas are the core
of the methodology presented in this section.

In the robust optimization framework, finding a safe approximation (or the
exact robust counterpart) for a problem that is convex (or concave) in the
uncertainty parameters relies on using convex (or concave) conjugate functions
to make the problem linear in the uncertainty (Bertsimas and den Hertog,
2022). The following lemmas are a generalization of this approach for the case
in which the objective involves the composition of two functions, with the
goal of making the problem linear not in the uncertainty but in the inner most
function. Lemma 1 makes this generalization when the outer function is convex
while Lemma 2 focuses on the concave case.

###### Lemma 1

If $f:A\rightarrow B$ is a convex and closed function then for any function
${\boldsymbol{z}}:\mathcal{U}\rightarrow{A}$, and any function $g:A\rightarrow
B$ we have

$\sup_{\boldsymbol{\delta}\in\mathcal{U}}\>f(\boldsymbol{z}(\boldsymbol{\delta}))+g(\boldsymbol{z}(\boldsymbol{\delta}))=\sup_{{\boldsymbol{u}}\in\text{dom}(f^{\star})}\enspace\sup_{\boldsymbol{\delta}\in\mathcal{U}}\>\boldsymbol{z}(\boldsymbol{\delta})^{T}{\boldsymbol{u}}-f^{\star}({\boldsymbol{u}})+g(\boldsymbol{z}(\boldsymbol{\delta})),$

where the convex conjugate function $f^{\star}$ is defined by
$f^{\star}({\boldsymbol{z}})=\sup_{{\boldsymbol{x}}\in\text{dom}(f)}\>{\boldsymbol{z}}^{T}{\boldsymbol{x}}-f({\boldsymbol{x}})$.

Proof  See Appendix A1.

###### Lemma 2

Let $g:A\rightarrow B$ be a concave and closed function. If a function
${\boldsymbol{z}}:\mathcal{U}\rightarrow{A}$ satisfies
$g({\boldsymbol{z}}(\boldsymbol{\delta}))<\infty$ for all
$\boldsymbol{\delta}\in\mathcal{U}$, then

$\sup_{\boldsymbol{\delta}\in\mathcal{U}}\>g({\boldsymbol{z}}(\boldsymbol{\delta}))=\inf_{{\boldsymbol{v}}\in\text{dom}(g_{\star})}\enspace\sup_{\boldsymbol{\delta}\in\mathcal{U}}\>{\boldsymbol{z}}(\boldsymbol{\delta})^{T}{\boldsymbol{v}}-g_{\star}({\boldsymbol{v}}),$

where the concave conjugate function is defined by
$g_{\star}({\boldsymbol{z}})=\inf_{{\boldsymbol{x}}\in\text{dom}(g)}\>{\boldsymbol{z}}^{T}{\boldsymbol{x}}-g({\boldsymbol{x}})$.

Proof  See Appendix A2.

From the lemmas above we can observe that to apply them we will need to
compute convex and concave conjugate functions. The next lemma facilitates
this computations for neural networks with ReLU activation functions.

###### Lemma 3

If ${\boldsymbol{u}},{\boldsymbol{p}},{\boldsymbol{q}}\geq\boldsymbol{0}$,
then the functions
$f({\boldsymbol{x}})={\boldsymbol{p}}^{\top}[{\boldsymbol{x}}]^{+}$ and
$g({\boldsymbol{x}})={\boldsymbol{x}}^{T}{\boldsymbol{u}}-{\boldsymbol{q}}^{\top}[{\boldsymbol{x}}]^{+}$
satisfy

$\displaystyle a)\>f^{\star}({\boldsymbol{z}})=\begin{cases}0&\quad\text{if
}\boldsymbol{0}\leq{\boldsymbol{z}}\leq{\boldsymbol{p}},\\\
\infty&\quad\text{otherwise.}\end{cases}\quad\text{and}\quad
b)\>g_{\star}({\boldsymbol{z}})=\begin{cases}0&\quad\text{if
}{\boldsymbol{u}}-{\boldsymbol{q}}\leq{\boldsymbol{z}}\leq{\boldsymbol{u}},\\\
-\infty&\quad\text{otherwise.}\end{cases}$

Proof  See Appendix A3.

As observed in Eqs. (7) to (9), we can obtain an upper bound of the min-max
problem in Eq. (2) by finding instead upper bounds for
$\max_{\boldsymbol{\delta}\in\mathcal{U}}\>(\boldsymbol{e}_{k}-{\boldsymbol{e}}_{y})^{\top}\boldsymbol{z}^{L}(\theta,{\boldsymbol{x}}+\boldsymbol{\delta})$
for each class $k$. We will find these upper bounds by recursively applying
the previous lemmas in each layer of the network. Starting from the last
layer, at each stage we split the model as the sum of a convex function and a
concave function, and we then apply Lemmas 1 and 2 to move the supremum over
$\mathcal{U}$ towards the previous layer until we reach the first layer of the
network. Since this first layer is linear in the uncertainty and the dual
function of the $L_{1}$ norm is the $L_{\infty}$ norm (a maximum over a finite
set), we are able to exactly solve the maximization problem over $\mathcal{U}$
and to backtrack to solve for the variables
${\boldsymbol{u}},{\boldsymbol{v}}$ that lemmas 1 and 2 introduce.

For simplicity, we develop and prove the upper bound for the robust
counterpart assuming that the neural network has only two layers; however the
results can be extended to the general case as shown in Appendix B. In
addition, all the theorems can be generalized for convolutions neural
networks, since convolutions are a special case of matrix multiplication.

The following theorem, whose proof relies on Lemmas 1 and 2, shows how the
supremum over $\mathcal{U}$ can be moved from the second layer to the first
layer.

###### Theorem 4

The maximum difference between the output of the correct class and the output
of any other class $k$ can be written as

$\displaystyle\quad\sup_{{\boldsymbol{{\delta}}}\in\mathcal{U}}\>{\boldsymbol{c}}_{k}^{\top}{\boldsymbol{z}}^{2}(\theta,{\boldsymbol{x}}+\boldsymbol{{\delta}})$
(12) $\displaystyle\begin{split}=&\sup_{0\leq{{\boldsymbol{s}}}\leq
1}\inf_{0\leq{{\boldsymbol{t}}}\leq
1}\sup_{\boldsymbol{{\delta}}\in\mathcal{U}}\>({\boldsymbol{p}}-{\boldsymbol{q}})^{\top}\>{\boldsymbol{z}}^{1}(\theta,{\boldsymbol{x}}+{\boldsymbol{{\delta}}})+{\boldsymbol{c}}_{k}^{\top}{\boldsymbol{b}}^{2}\\\
&\hskip
42.67912pt\text{s.t.}\quad{\boldsymbol{p}}=[({\boldsymbol{W}}^{2})^{\top}{\boldsymbol{c}}_{k}]^{+}\odot{\boldsymbol{s}}\\\
&\hskip
68.28644pt{\boldsymbol{q}}=[-({\boldsymbol{W}}^{2})^{\top}{\boldsymbol{c}}_{k}]^{+}\odot{\boldsymbol{t}},\\\
\end{split}$ (13)

where $\odot$ corresponds to entry-wise multiplication.

Proof  By definition of the two layer neural network ${\boldsymbol{z}}^{2}$,
we have

$\displaystyle{\boldsymbol{c}}_{k}^{\top}{\boldsymbol{z}}^{2}(\theta,{\boldsymbol{x}}+{\boldsymbol{\delta}})$
$\displaystyle={\boldsymbol{c}}_{k}^{\top}{\boldsymbol{W}}^{2}[{\boldsymbol{z}}^{1}(\theta,{\boldsymbol{x}}+{\boldsymbol{\delta}})]^{+}+{\boldsymbol{c}}_{k}^{\top}{\boldsymbol{b}}^{2}$
$\displaystyle=f_{+}({\boldsymbol{z}}^{1}(\theta,{\boldsymbol{x}}+{\boldsymbol{\delta}}))-f_{-}({\boldsymbol{z}}^{1}(\theta,{\boldsymbol{x}}+{\boldsymbol{\delta}}))+{\boldsymbol{c}}_{k}^{\top}{\boldsymbol{b}}^{2},$

where $f_{+},f_{-}$ are the convex functions defined by

$f_{+}({\boldsymbol{x}})=[{\boldsymbol{c}}_{k}^{\top}{\boldsymbol{W}}^{2}]^{+}[{\boldsymbol{x}}]^{+},\quad\text{
and }\quad
f_{-}({\boldsymbol{x}})=[-{\boldsymbol{c}}_{k}^{\top}{\boldsymbol{W}}^{2}]^{+}[{\boldsymbol{x}}]^{+}.$

Applying Lemma 1 to the function $f_{+}$ we then have

$\displaystyle\sup_{{\boldsymbol{\delta}}\in\mathcal{U}}\>{\boldsymbol{c}}_{k}^{\top}{\boldsymbol{z}}^{2}(\theta,{\boldsymbol{x}}+{\boldsymbol{\delta}}),$
$\displaystyle=$
$\displaystyle\sup_{{\boldsymbol{\delta}}\in\mathcal{U}}\>f_{+}({\boldsymbol{z}}^{1}(\theta,{\boldsymbol{x}}+{\boldsymbol{\delta}}))-f_{-}({\boldsymbol{z}}^{1}(\theta,{\boldsymbol{x}}+{\boldsymbol{\delta}}))+{\boldsymbol{c}}_{k}^{\top}{\boldsymbol{b}}^{2},$
$\displaystyle=$
$\displaystyle\sup_{{\boldsymbol{u}}\in\text{dom}(f_{+}^{\star})}\>\sup_{{\boldsymbol{\delta}}\in\mathcal{U}}\>\>{\boldsymbol{u}}^{\top}{\boldsymbol{z}}^{1}(\theta,{\boldsymbol{x}}+{\boldsymbol{\delta}})-f_{+}^{\star}({\boldsymbol{u}})-f_{-}({\boldsymbol{z}}^{1}(\theta,{\boldsymbol{x}}+{\boldsymbol{\delta}}))+{\boldsymbol{c}}_{k}^{\top}{\boldsymbol{b}}^{2},\quad$
$\displaystyle(\text{By Lemma \ref{convex-bound}}),$ (14) $\displaystyle=$
$\displaystyle\sup_{{\boldsymbol{u}}\in\text{dom}(f_{+}^{\star})}\>\sup_{{\boldsymbol{\delta}}\in\mathcal{U}}\>\>{\boldsymbol{u}}^{\top}{\boldsymbol{z}}^{1}(\theta,{\boldsymbol{x}}+{\boldsymbol{\delta}})-f_{-}({\boldsymbol{z}}^{1}(\theta,{\boldsymbol{x}}+{\boldsymbol{\delta}}))+{\boldsymbol{c}}_{k}^{\top}{\boldsymbol{b}}^{2},\quad$
$\displaystyle(\text{By Lemma \ref{conjugate-computation}a)}.$ (15)

Defining the concave function
$g({\boldsymbol{x}})={\boldsymbol{u}}^{\top}{\boldsymbol{x}}-f_{-}({\boldsymbol{x}})$,
and applying Lemma 2 to the function $g$ we obtain

$\displaystyle\sup_{{\boldsymbol{\delta}}\in\mathcal{U}}\>{\boldsymbol{c}}_{k}^{\top}{\boldsymbol{z}}^{2}(\theta,{\boldsymbol{x}}+{\boldsymbol{\delta}})$
$\displaystyle=$
$\displaystyle\sup_{{\boldsymbol{u}}\in\text{dom}(f_{+}^{\star})}\>\sup_{{\boldsymbol{\delta}}\in\mathcal{U}}\>g({\boldsymbol{z}}^{1}(\theta,{\boldsymbol{x}}+\boldsymbol{\delta}))+{\boldsymbol{c}}_{k}^{\top}{\boldsymbol{b}}^{2}$
(16) $\displaystyle=$
$\displaystyle\sup_{{\boldsymbol{u}}\in\text{dom}(f_{+}^{\star})}\>\inf_{{\boldsymbol{v}}\in\text{dom}(g_{\star})}\>\sup_{{\boldsymbol{\delta}}\in\mathcal{U}}\>{\boldsymbol{v}}^{\top}{\boldsymbol{z}}^{1}(\theta,{\boldsymbol{x}}+{\boldsymbol{\delta}})-g_{\star}({\boldsymbol{v}})+{\boldsymbol{c}}_{k}^{\top}{\boldsymbol{b}}^{2}\quad$
$\displaystyle(\text{By Lemma \ref{concave-bound}}),$ (17) $\displaystyle=$
$\displaystyle\sup_{{\boldsymbol{u}}\in\text{dom}(f_{+}^{\star})}\>\inf_{{\boldsymbol{v}}\in\text{dom}(g_{\star})}\>\sup_{{\boldsymbol{\delta}}\in\mathcal{U}}\>{\boldsymbol{v}}^{\top}{\boldsymbol{z}}^{1}(\theta,{\boldsymbol{x}}+{\boldsymbol{\delta}})+{\boldsymbol{c}}_{k}^{\top}{\boldsymbol{b}}^{2}\quad$
$\displaystyle(\text{By Lemma \ref{conjugate-computation}b}).$ (18)

Lastly, by Lemma 3a and 3b, we know that the variables ${\boldsymbol{u}}$ and
${\boldsymbol{v}}$ can be parameterized as

$\displaystyle{\boldsymbol{u}}$
$\displaystyle=[({\boldsymbol{W}}^{2})^{\top}{\boldsymbol{c}}_{k}]^{+}\odot{\boldsymbol{s}},$
$\displaystyle{\boldsymbol{v}}$
$\displaystyle=[({\boldsymbol{W}}^{2})^{\top}{\boldsymbol{c}}_{k}]^{+}\odot{\boldsymbol{s}}-[-({\boldsymbol{W}}^{2})^{\top}{\boldsymbol{c}}_{k}]^{+}\odot{\boldsymbol{t}}$

with $0\leq{\boldsymbol{s}},{\boldsymbol{t}}\leq 1$. Substituting these values
in Eq. (18) we obtain Eq. (LABEL:induction-problem1), as desired.

Notice that the objective in Eq. (LABEL:induction-problem1) is linear in
${\boldsymbol{z}}^{1}$ and therefore it is also linear in
${\boldsymbol{\delta}}$, which facilitates the computation of the exact value
of the supremum over $\mathcal{U}$, as shown in the next corollary.

###### Corollary 5

If
$\mathcal{U}=\\{\boldsymbol{\delta}:\|\boldsymbol{\delta}\|_{p}\leq\rho\\}$,
then:

$\displaystyle\sup_{{\boldsymbol{\delta}}\in\mathcal{U}}\>{\boldsymbol{c}}_{k}^{\top}{\boldsymbol{z}}^{2}(\theta,{\boldsymbol{x}}+{\boldsymbol{\delta}})$
(19) $\displaystyle\begin{split}=&\sup_{0\leq{{\boldsymbol{s}}}\leq
1}\inf_{0\leq{{\boldsymbol{t}}}\leq
1}\>\rho\|({\boldsymbol{p}}-{\boldsymbol{q}})^{\top}{\boldsymbol{W}}^{1}\|_{q}+({\boldsymbol{p}}-{\boldsymbol{q}})^{\top}({\boldsymbol{W}}^{1}{\boldsymbol{x}}+{\boldsymbol{b}}^{1})+{\boldsymbol{c}}_{k}^{\top}{\boldsymbol{b}}^{2}\\\
&\hskip
42.67912pt\text{s.t.}\quad{\boldsymbol{p}}=[({\boldsymbol{W}}^{2})^{\top}{\boldsymbol{c}}_{k}]^{+}\odot{\boldsymbol{s}}\\\
&\hskip
68.28644pt{\boldsymbol{q}}=[-({\boldsymbol{W}}^{2})^{\top}{\boldsymbol{c}}_{k}]^{+}\odot{\boldsymbol{t}},\\\
\end{split}$ (20)

where $\|\cdot\|_{q}$ is the dual norm of $\|\cdot\|_{p}$, with
$\frac{1}{p}+\frac{1}{q}=1$.

Before proceeding to the proof of the corollary, notice that we can recover
the approximation method developed in the previous section by setting

$\displaystyle{\boldsymbol{s}}={\boldsymbol{t}}=[sign({\boldsymbol{z}}^{1}(\theta,{\boldsymbol{x}}))]^{+}$
(21)

in the objective of problem (LABEL:conjugate_problem1) to obtain

$\displaystyle\rho\|(({\boldsymbol{W}}^{2})^{\top}{\boldsymbol{c}}_{k}\odot[sign({\boldsymbol{z}}^{1}(\theta,{\boldsymbol{x}}))]^{+})^{\top}{\boldsymbol{W}}^{1}\|_{q}+{\boldsymbol{c}}_{k}^{\top}{\boldsymbol{W}}^{2}[{\boldsymbol{z}}^{1}(\theta,{\boldsymbol{x}})]^{+}+{\boldsymbol{c}}_{k}^{\top}{\boldsymbol{b}}^{2},$
(22)

which is the same as the linear approximation of
${\boldsymbol{c}}_{k}^{\top}{\boldsymbol{z}}^{2}(\theta+\boldsymbol{\delta})$
obtained in Eq. (5).

Proof  The proof follows directly after applying Theorem 4 and using the fact
that for all vectors $\boldsymbol{{\boldsymbol{c}}}$ we have

$\displaystyle\sup_{{\boldsymbol{\delta}}:\|\boldsymbol{\delta}\|_{p}\leq\rho}{\boldsymbol{c}}^{\top}\boldsymbol{\delta}=\rho\|{\boldsymbol{c}}\|_{q}.$
(23)

Since neural networks are trained by minimizing the empirical loss over the
parameters $\theta$, we want to avoid the computation of supremums in the
objective. While the previous corollary shows how to solve the supremum over
the uncertainty set, a new supremum was introduced in Theorem 4 over the
variables ${\boldsymbol{s}}$. The next theorem tells us how we can remove this
new supremums for the specific case $p=1$.

###### Theorem 6

The maximum difference between the output of the correct class and the output
of any other class $k$ can be upper bounded by

$\displaystyle\sup_{{\boldsymbol{\delta}}:\|\boldsymbol{\delta}\|_{1}\leq\rho}$
$\displaystyle\>{\boldsymbol{c}}_{k}^{\top}{\boldsymbol{z}}^{2}(\theta,{\boldsymbol{x}}+{\boldsymbol{\delta}})$
$\displaystyle\leq$ $\displaystyle\inf_{0\leq{{\boldsymbol{t}}}\leq
1}\>\max_{m\in[M]}\>\max\bigg{\\{}\
g^{2}_{k,m}(\theta,{\boldsymbol{t}},{\boldsymbol{x}},\rho),g^{2}_{k,m}(\theta,{\boldsymbol{t}},{\boldsymbol{x}},-\rho)\bigg{\\}},$
(24)

where the new network $g$ is defined by the equations

$\displaystyle g^{1}_{m}(\theta,a)$
$\displaystyle=a{\boldsymbol{W}}^{1}\boldsymbol{e}_{m}+{\boldsymbol{W}}^{1}{\boldsymbol{x}}+{\boldsymbol{b}}^{1},$
$\displaystyle g^{2}_{k,m}(\theta,{\boldsymbol{t}},{\boldsymbol{x}},a)$
$\displaystyle=[{\boldsymbol{c}}_{k}^{\top}{\boldsymbol{W}}^{2}]^{+}[g^{1}_{m}(\theta,a)]^{+}-[-{\boldsymbol{c}}_{k}^{\top}{\boldsymbol{W}}^{2}]^{+}[g^{1}_{m}(\theta,a)]\odot{\boldsymbol{t}}+{\boldsymbol{c}}_{k}^{\top}{\boldsymbol{b}}^{2},$

for $a=\rho,-\rho$.

Proof  Applying Corollary 5 with $p=1$ and using the min-max inequality we
obtain

$\displaystyle\sup_{{\boldsymbol{\delta}}:\|{\boldsymbol{\delta}}\|_{1}\leq\rho}\>{\boldsymbol{c}}_{k}^{\top}{\boldsymbol{z}}^{2}(\theta,{\boldsymbol{x}}+{\boldsymbol{\delta}})$
(25) $\displaystyle\begin{split}\leq&\inf_{0\leq{{\boldsymbol{t}}}\leq
1}\sup_{0\leq{{\boldsymbol{s}}}\leq
1}\rho\|({\boldsymbol{p}}-{\boldsymbol{q}})^{\top}{\boldsymbol{W}}^{1}\|_{\infty}+({\boldsymbol{p}}-{\boldsymbol{q}})^{\top}({\boldsymbol{W}}^{1}{\boldsymbol{x}}+{\boldsymbol{b}}^{1})+{\boldsymbol{c}}_{k}^{\top}{\boldsymbol{b}}^{2}\\\
&\hskip
42.67912pt\text{s.t.}\quad{\boldsymbol{p}}=[({\boldsymbol{W}}^{2})^{\top}{\boldsymbol{c}}_{k}]^{+}\odot{\boldsymbol{s}}\\\
&\hskip
68.28644pt{\boldsymbol{q}}=[-({\boldsymbol{W}}^{2})^{\top}{\boldsymbol{c}}_{k}]^{+}\odot{\boldsymbol{t}}.\end{split}$
(26)

Defining
${\boldsymbol{p}}({\boldsymbol{s}})=[({\boldsymbol{W}}^{2})^{\top}{\boldsymbol{c}}_{k}]^{+}\odot{\boldsymbol{s}}$
and
${\boldsymbol{q}}({\boldsymbol{t}})=[-({\boldsymbol{W}}^{2})^{\top}{\boldsymbol{c}}_{k}]^{+}\odot{\boldsymbol{t}}$,
we have that for fixed ${{\boldsymbol{t}}}$ it holds

$\displaystyle\sup_{0\leq{{\boldsymbol{s}}}\leq
1}\>\rho\|({\boldsymbol{p}}({\boldsymbol{s}})-{\boldsymbol{q}}({\boldsymbol{t}}))^{\top}{\boldsymbol{W}}^{1}\|_{\infty}+({\boldsymbol{p}}({\boldsymbol{s}})-{\boldsymbol{q}}({\boldsymbol{t}}))^{\top}({\boldsymbol{W}}^{1}{\boldsymbol{x}}+{\boldsymbol{b}}^{1})+{\boldsymbol{c}}_{k}^{\top}{\boldsymbol{b}}^{2}$
$\displaystyle=$
$\displaystyle\max_{m\in[M]}\>\max\bigg{\\{}\sup_{0\leq{{\boldsymbol{s}}}\leq
1}\>({\boldsymbol{p}}({\boldsymbol{s}})-{\boldsymbol{q}}({\boldsymbol{t}}))^{\top}({\boldsymbol{W}}^{1}({\boldsymbol{x}}+\rho{\boldsymbol{e}}_{m})+{\boldsymbol{b}}^{1})+{\boldsymbol{c}}_{k}^{\top}{\boldsymbol{b}}^{2},$
$\displaystyle\hskip 59.75095pt\sup_{0\leq{{\boldsymbol{s}}}\leq
1}\>({\boldsymbol{p}}({\boldsymbol{s}})-{\boldsymbol{q}}({\boldsymbol{t}}))^{\top}({\boldsymbol{W}}^{1}({\boldsymbol{x}}-\rho{\boldsymbol{e}}_{m})+{\boldsymbol{b}}^{1})+{\boldsymbol{c}}_{k}^{\top}{\boldsymbol{b}}^{2}\bigg{\\}}$
$\displaystyle=$ $\displaystyle\max_{m\in[M]}\>\max\bigg{\\{}\
[{\boldsymbol{c}}_{k}^{\top}{\boldsymbol{W}}^{2}]^{+}[{\boldsymbol{W}}^{1}({\boldsymbol{x}}+\rho{\boldsymbol{e}}_{m})+{\boldsymbol{b}}^{1}]^{+}-{\boldsymbol{q}}({\boldsymbol{t}})^{\top}({\boldsymbol{W}}^{1}({\boldsymbol{x}}+\rho{\boldsymbol{e}}_{m})+{\boldsymbol{b}}^{1})+{\boldsymbol{c}}_{k}^{\top}{\boldsymbol{b}}^{2},$
$\displaystyle\hskip
59.75095pt[({\boldsymbol{c}}_{k}^{\top}{\boldsymbol{W}}^{2}]^{+}[{\boldsymbol{W}}^{1}({\boldsymbol{x}}-\rho{\boldsymbol{e}}_{m})+{\boldsymbol{b}}^{1}]^{+}-{\boldsymbol{q}}({\boldsymbol{t}})^{\top}({\boldsymbol{W}}^{1}({\boldsymbol{x}}-\rho{\boldsymbol{e}}_{m})+{\boldsymbol{b}}^{1})+{\boldsymbol{c}}_{k}^{\top}{\boldsymbol{b}}^{2}\bigg{\\}}$
$\displaystyle=$
$\displaystyle\max_{m\in[M]}\>\max\\{g^{2}_{k,m}(\theta,{\boldsymbol{t}},{\boldsymbol{x}},\rho),g^{2}_{k,m}(\theta,{\boldsymbol{t}},{\boldsymbol{x}},-\rho)\\}.$

The theorem then follows after applying the $\inf$ over
$0\leq{\boldsymbol{t}}\leq 1$.
Notice that in the previous proof it was crucial to use $p=1$, since the dual
of the $L_{1}$ norm is the $L_{\infty}$ norm, which can be written as a
maximum over a finite set. With a different $p$ we would not have been able to
apply the same technique to remove the variables ${\boldsymbol{s}}$. However,
for the chosen uncertainty set we obtain an upper bound of Eq. (2) by applying
the result from the previous Theorem to Eq. (9).

While we could include the variables ${\boldsymbol{t}}$ in the minimization
problem over $\theta$, we instead use fixed values
${\boldsymbol{t}}=[sign({\boldsymbol{z}}^{1}(\theta,{\boldsymbol{x}}))]^{+}$
based on the linear approximation of
${\boldsymbol{c}}_{k}^{\top}{\boldsymbol{z}}^{2}(\theta,{\boldsymbol{x}}+{\boldsymbol{\delta}})$,
as described in Eq. (21). Notice that setting specific values for
${\boldsymbol{t}}$ does not affect the inequalities: since the upper bound
includes the infimum over ${\boldsymbol{t}}$, any $0\leq{\boldsymbol{t}}\leq
1$ yields an upper bound of the robust problem. For the specific case of the
cross entropy loss function, the proposed upper bound for the min-max robust
problem is

$\displaystyle\min_{\theta}\frac{1}{N}\sum_{n=1}^{N}\log\bigg{(}\sum_{k}e^{\big{(}{\max_{m\in[M]}\max\\{g^{2}_{k,m}(\theta,[sign({\boldsymbol{z}}^{1}(\theta,{\boldsymbol{x}}))]^{+},{\boldsymbol{x}},\rho),g^{2}_{k,m}(\theta,[sign({\boldsymbol{z}}^{1}(\theta,{\boldsymbol{x}}))]^{+},{\boldsymbol{x}},-\rho)\\}}\big{)}}\bigg{)}.$
(27)

##  6 Experiments

In this section, we demonstrate the effectiveness of the proposed methods in
practice. We first introduce the experimental setup and then we compare the
robustness of several defenses.

###  6.1 Experimental details

#### Data sets.

We use 46 data sets from the UCI collection (Dua and Graff, 2017), which
correspond to classification tasks with a diverse number of features that are
not categorical. For each data set we do a $80\%/20\%$ split for
training/testing sets, and we further reserve $25\%$ of each training set for
validation. In addition, we use three popular computer vision data sets,
namely the MNIST (Deng, 2012), Fashion MNIST (Xiao et al., 2017) and CIFAR
(Krizhevsky, 2009) data sets.

#### Pre-processing.

All input data has been previously scaled, which facilitates the comparison of
the adversarial attacks across data sets. For the UCI data sets, each feature
is standarized using the statistics of the training set, while for the vision
data sets each image channel is normalized to be between $0$ and $1$, or
standarized, depending on what leads to best robustness.

#### Attacks.

We use the implementation provided by the foolbox library (Rauber et al.,
2017, 2020) using the default parameters. We evaluate the PGD, FGM and FGSM
attacks. More specifically, we use the following attacks:

  * •

PGD-$L_{p}$: Adversarial attack bounded in $L_{p}$ norm and found using
Projected Gradient Descent.

  * •

FGM-$L_{p}$: Adversarial attack bounded in $L_{p}$ norm and found using Fast
Gradient Method.

  * •

FSGM-$L_{p}$: Adversarial attack bounded in $L_{p}$ norm and found using Fast
Sign Gradient Method.

For the UCI data sets we evaluate adversarial accuracies with respect to
attacks bounded in $L_{1},L_{2}$ and $L_{\infty}$ norm. For the vision data
sets, we present results only for attacks bounded in the $L_{2}$ norm, but
results using other attacks can be found in Appendix C.

#### Defenses.

Our comparisons include different defenses denoted as follows:

  * •

ARUB-$L_{p}$: Approximate Robust Upper Bound method described in section 4
using the $L_{1}$ or $L_{\infty}$ sphere as the uncertainty set.

  * •

RUB: Robust Upper Bound method described in section 5 using the $L_{1}$ sphere
as the uncertainty set.

  * •

PGD-$L_{\infty}$: Adversarial training method in which the attacks are bounded
in the $L_{\infty}$ norm and found using Projected Gradient Descent.

  * •

Nominal: Standard vanilla training with no robustness ($\rho=0$).

  * •

Baseline-$L_{\infty}$: Simple approximation method resulting from minimizing
Eq. (3) using the $L_{p}$ sphere as the uncertainty set.

#### Architecture.

We evaluate a neural network with three dense hidden layers with 200 neurons
in each hidden layer. For the vision datasets, we also provide results with
Convolutional Neural Networks (_CNNs_) in Appendix D. The architecture has two
convolutional layers alternated with pooling operations, and two dense layers,
as in Madry et al. (2019).

#### Hyperparameter Tuning.

Each network and defense is trained multiple times for different learning
rates ($\\{1,10^{-1},10^{-2},10^{-3},10^{-4},10^{-5},10^{-6}\\}$). For the UCI
data set we use a batch size of $256$ and for the vision data sets we try a
batch size of $32$ and $256$. All networks trained using the UCI data sets are
trained for $5000$ iterations, and all vision data sets are trained for
$10000$ iteration which is sufficient for convergence in both cases. For the
$L_{\infty}$ based training methods we try all values of $\rho$ from the set
($\\{10^{-4},10^{-5},10^{-3},10^{-2},0.1,0.3$, $0.5,1,$ $3,5,10\\}$). For the
methods based on the $L_{1}$ norm, we scale those values of $\rho$ by a factor
of $\sqrt{m}$, since
$\|{\boldsymbol{x}}\|_{\infty}\leq\|{\boldsymbol{x}}\|_{2}$ and
$\|{\boldsymbol{x}}\|_{1}\leq\sqrt{m}\|{\boldsymbol{x}}\|_{2}$ for any
${\boldsymbol{x}}\in\mathbb{R}^{m}$. In this way we ensure that the $L_{1}$
spheres and the $L_{\infty}$ spheres contain the same $L_{2}$ spheres,
allowing for a fair comparison of all methods in terms of adversarial attacks
that are bounded in the $L_{2}$ norm. For each attack and each $\rho$ of the
attack, we select on the validation set the best hyperparameters for each
defense.

###  6.2 UCI data sets

We run experiments on the $46$ UCI data sets using different methods for
robust training and compare the adversarial accuracies achieved with multiple
types of adversarial attacks. For each data set we rank every training method,
where the method with rank $1$ corresponds to the one with highest adversarial
accuracy. The results are shown in Figure 2, where we observe a similar
pattern across all types of attacks, namely, we see that the best ranks are
achieved with ARUB-$L_{\infty}$ and PGD-$L_{\infty}$ when $\rho$ is smaller
than $10^{-1}$; next there is a small range in which PGD-$L_{\infty}$ does
best and finally for larger values of $\rho$ the best rank is that of RUB. We
also highlight that looking at $\rho=0$, it is clear that robust training
methods achieve better natural accuracy than the one obtained with Nominal
training.

To better compare the performances of RUB and PGD-$L_{\infty}$, in Figure 3 we
show the number of data sets for which RUB improves $L_{2}$ adversarial
accuracy over PGD-$L_{\infty}$ by a specific percentage. We observe that the
improvement becomes more significant as $\rho$ increases, and in particular,
for $\rho=10$ we observe that RUB only lowers adversarial accuracy for $3$
data sets, while it shows more than $15\%$ improvement over PGD-$L_{\infty}$
for $8$ of the data sets.

(a)

(b)

(c)

(d)

(e)

(f)

Figure 2: Average rank of each method accross the 46 UCI data sets for
adversarial attacks bounded in $L_{2},L_{\infty}$ and $L_{1}$ norm,
respectively from top to bottom. The figures on the left use attacks based on
Fast Gradient methods, and the figures on the right use instead attacks based
on Projected Gradient Descent.

(a) $\rho=0.01$

(b) $\rho=1$

(c) $\rho=10$

Figure 3: Number of UCI data sets for which RUB-$L_{1}$ improves adversarial
accuracy over PGD-$L_{\infty}$ by a specific percentage. Figures a), b) and c)
show the corresponding plots for PGD-$L_{2}$ adversarial examples with
different values of $\rho$.

Lastly, in Table 2 we display the total training time of each method in a
$\log$ scale. As expected, we see that nominal is the fastest training method,
and while ARUB-$L_{\infty}$ and RUB increase the total time by several orders
of magnitude, they are significantly faster than PGD-$L_{\infty}$.

| Avg no. samples per second | Standard Deviation
---|---|---
RUB | 9.71 | 8.07
ARUB -$L_{\infty}$ | 8.47 | 4.74
Baseline-$L_{\infty}$ | 11.69 | 9.58
PGD -$L_{\infty}$ | 7.06 | 2.60
Nominal | 12.11 | 10.02
Table 2: Average $\log$ of the number of samples processed per second across
the 46 UCI data sets, as well as the corresponding $\log$ of the standard
deviations.

###  6.3 Vision Data Sets

We next show experiment results for the three vision data sets. Specifically,
we compare for different training methods their performance against
adversarial attacks as well as the security guarantees obtained from applying
the upper bound from Eq. (27). We report adversarial accuracies using
PGD-$L_{2}$ attacks. Also, adversarial accuracies using using attacks bounded
in $L_{1}$ and $L_{\infty}$ norms can be found in Appendix C. Moreover, since
the proposed RUB method significantly increases memory requirements for inputs
with large dimensions, we compare this method against other defenses using
only the feed architecture with three hidden layers. We do include results
using a CNN architecture for all other methods in Appendix D, where we also
explain how to extend the theory of the RUB method for convolutional layers
with ReLU and MaxPool activation functions.

#### Performance Against Adversarial Attacks.

We evaluate adversarial accuracy for 6 different methods: RUB, ARUB-$L_{1}$,
ARUB-$L_{\infty}$, Baseline-$L_{\infty}$, PGD-$L_{\infty}$, and Nominal
training. We observe that for the Fashion MNIST data set (Table 3) as well as
the MNIST data set (Table 4), ARUB-$L_{1}$ achieves the best accuracies for
small values of $\rho$; we then observe a small range in which
PGD-$L_{\infty}$ does best and lastly for larger values of $\rho$ we see that
RUB takes the lead, which is similar to the average results observed for the
UCI data sets. For the CIFAR data set we observe a different behavior;
PGD-$L_{\infty}$ achieves better accuracies for small and large values of
$\rho$, although we observe that all methods perform very poorly overall as
this is a notoriously more difficult data set. Lastly, once again we observe
that in all three data sets the robust training methods achieve a higher
natural accuracy than the accuracy resulting from standard nominal training.

$\boldsymbol{\rho}=$ | 0.00 | 0.01 | 0.06 | 0.28 | 2.80 | 8.40 | 14.00 | 28.00
---|---|---|---|---|---|---|---|---
RUB | 90.59 | 89.88 | 88.95 | 86.80 | 63.09 | 56.95 | 57.07 | 55.27
ARUB-$L_{1}$ | 90.31 | 90.04 | 89.45 | 85.98 | 66.17 | 30.63 | 18.24 | 15.04
ARUB-$L_{\infty}$ | 89.18 | 89.06 | 88.75 | 85.98 | 67.11 | 24.18 | 18.24 | 16.02
Baseline-$L_{\infty}$ | 89.38 | 89.02 | 88.32 | 85.39 | 50.20 | 29.88 | 28.71 | 26.52
PGD-$L_{\infty}$ | 89.61 | 89.38 | 88.01 | 85.23 | 68.20 | 31.60 | 25.90 | 22.85
Nominal | 87.70 | 88.44 | 88.01 | 85.35 | 46.60 | 19.92 | 16.84 | 15.39
Table 3: Adversarial Accuracy for Fashion MNIST with PGD-$L_{2}$ attacks $\boldsymbol{\rho}=$ | 0.00 | 0.01 | 0.06 | 0.28 | 2.80 | 8.40 | 14.00 | 28.00
---|---|---|---|---|---|---|---|---
RUB | 97.93 | 97.89 | 97.07 | 97.50 | 79.80 | 66.48 | 62.89 | 55.00
ARUB-$L_{1}$ | 98.63 | 98.52 | 97.73 | 97.46 | 81.13 | 57.58 | 55.04 | 46.60
ARUB-$L_{\infty}$ | 97.97 | 98.09 | 97.81 | 98.28 | 92.62 | 51.80 | 49.88 | 43.40
Baseline-$L_{\infty}$ | 97.58 | 97.50 | 97.58 | 97.11 | 72.85 | 64.18 | 61.29 | 52.81
PGD-$L_{\infty}$ | 98.24 | 98.32 | 98.24 | 98.09 | 93.12 | 66.05 | 62.46 | 54.26
Nominal | 97.62 | 97.46 | 97.34 | 96.52 | 69.80 | 19.57 | 18.63 | 17.70
Table 4: Adversarial Accuracy for MNIST with PGD-$L_{2}$ attacks $\boldsymbol{\rho}=$ | 0.00 | 0.01 | 0.06 | 0.55 | 5.54 | 16.63 | 27.71 | 55.43
---|---|---|---|---|---|---|---|---
RUB | 51.17 | 50.23 | 48.59 | 44.14 | 20.62 | 9.96 | 9.96 | 9.96
ARUB-$L_{1}$ | 52.70 | 51.41 | 51.13 | 46.99 | 23.79 | 9.69 | 9.30 | 9.30
ARUB-$L_{\infty}$ | 53.83 | 52.89 | 47.77 | 41.37 | 22.66 | 11.48 | 11.52 | 11.56
Baseline-$L_{\infty}$ | 53.12 | 52.07 | 46.52 | 37.19 | 9.30 | 9.30 | 9.30 | 9.30
PGD-$L_{\infty}$ | 53.91 | 53.32 | 48.83 | 45.66 | 22.19 | 13.63 | 13.48 | 11.99
Nominal | 46.99 | 46.88 | 45.59 | 38.59 | 10.59 | 10.59 | 10.59 | 10.59
Table 5: Adversarial Accuracy for CIFAR with PGD-$L_{2}$ attacks

#### Security Guarantees against $L_{1}$ Norm Bounded Attacks.

Finally, we use the upper bound of the robust loss derived in section 4 to
find lower bounds for the adversarial accuracy with respect to attacks bounded
in the $L_{1}$ norm by $\rho$. In Figures 8, 6 and 7 we observe that, as
expected, for all three data sets, CIFAR, Fashion MNIST and MNIST, the best
security guarantees are the ones for the RUB method. While these results are
only lower bounds for the adversarial accuracy and we cannot claim a better
accuracy for RUB than for the rest of the methods, the lower bound for the RUB
method shows that this method indeed performs very well against $L_{1}$
attacks bounded by large values of $\rho$. For instance, in Table 6 we can see
that for the Fashion MNIST data set the RUB method guarantees $86.02\%$
adversarial accuracy (less than $5\%$ decrease from the best natural accuracy)
against attacks with $L_{1}$ norm smaller or equal to $\rho=2.8$. Similarly,
in Table 7 we observe that for this same attacks RUB has at least $97.11\%$
adversarial accuracy (less than $1\%$ decrease over natural accuracy) for the
MNIST data set. And finally, for the CIFAR data set, we can see in Table 8
that RUB achieves $45.66\%$ adversarial accuracy (less than $5\%$ decrease
over natural accuracy) against attacks whose $L_{1}$ norm is upper bounded by
$\rho=5.54$.

$\boldsymbol{\rho}=$ | 0.00 | 0.01 | 0.06 | 0.28 | 2.80 | 8.40 | 14.00 | 28.00
---|---|---|---|---|---|---|---|---
RUB | 90.27 | 90.16 | 90.35 | 89.49 | 86.02 | 80.55 | 76.21 | 66.95
ARUB-$L_{1}$ | 90.51 | 90.51 | 90.39 | 89.73 | 85.59 | 73.98 | 69.61 | 47.54
ARUB-$L_{\infty}$ | 89.96 | 89.92 | 89.84 | 88.75 | 76.60 | 21.37 | 10.16 | 9.84
Baseline-$L_{\infty}$ | 89.49 | 89.38 | 90.16 | 87.50 | 77.42 | 46.41 | 20.04 | 15.23
PGD-$L_{\infty}$ | 89.92 | 90.74 | 90.47 | 87.85 | 78.75 | 40.35 | 19.22 | 15.55
Nominal | 88.59 | 88.55 | 88.48 | 87.93 | 77.50 | 40.98 | 15.04 | 9.88
Table 6: Lower Bound for the adversarial accuracy achieved for the Fashion MNIST data set with respect to uncertainty sets bounded in $L_{1}$ norm by $\rho$. $\boldsymbol{\rho}=$ | 0.00 | 0.01 | 0.06 | 0.28 | 2.80 | 8.40 | 14.00 | 28.00
---|---|---|---|---|---|---|---|---
RUB | 98.01 | 97.93 | 97.93 | 97.46 | 97.11 | 93.67 | 89.77 | 74.96
ARUB-$L_{1}$ | 98.52 | 98.48 | 98.48 | 98.20 | 96.48 | 89.96 | 80.86 | 26.05
ARUB-$L_{\infty}$ | 98.40 | 98.40 | 98.28 | 97.54 | 94.69 | 68.52 | 16.56 | 12.19
Baseline-$L_{\infty}$ | 98.05 | 98.05 | 98.05 | 97.81 | 95.23 | 57.23 | 20.74 | 10.35
PGD-$L_{\infty}$ | 98.48 | 98.48 | 98.44 | 98.36 | 95.55 | 63.59 | 15.43 | 11.60
Nominal | 97.73 | 97.73 | 97.58 | 97.54 | 93.87 | 44.80 | 10.31 | 10.23
Table 7: Lower Bound for the adversarial accuracy achieved for the MNIST data set with respect to uncertainty sets bounded in $L_{1}$ norm by $\rho$. $\boldsymbol{\rho}=$ | 0.00 | 0.01 | 0.06 | 0.55 | 5.54 | 16.63 | 27.71 | 55.43
---|---|---|---|---|---|---|---|---
RUB | 50.62 | 50.62 | 49.96 | 48.91 | 45.66 | 37.81 | 32.85 | 23.28
ARUB-$L_{1}$ | 53.40 | 53.16 | 51.99 | 51.33 | 43.36 | 33.83 | 27.19 | 15.16
ARUB-$L_{\infty}$ | 53.67 | 53.52 | 53.20 | 47.07 | 37.30 | 14.26 | 9.65 | 9.65
Baseline-$L_{\infty}$ | 53.32 | 53.24 | 52.66 | 46.09 | 36.99 | 13.71 | 9.61 | 9.61
PGD-$L_{\infty}$ | 54.88 | 54.80 | 53.98 | 49.26 | 37.42 | 27.30 | 14.02 | 12.46
Nominal | 46.88 | 46.88 | 46.76 | 44.69 | 37.19 | 14.88 | 9.02 | 8.95
Table 8: Lower Bound for the adversarial accuracy achieved for the CIFAR data
set with respect to uncertainty sets bounded in $L_{1}$ norm by $\rho$.

##  7 Conclusions

We developed two new methods for robust training of neural networks, both of
which provide a new objective function in closed form solution that can be
effectively trained with backpropagation. First, we found an approximated
upper bound for the robust loss by incorporating the first order approximation
of the network’s output layer. This method does not provide security
guarantees against adversarial attacks but it performs very well across a
variety of data sets when the uncertainty set is small. Second, by extending
state-of-the-art tools from robust optimization to non-convex and non-concave
functions, we were able to construct an exact upper bound for the robust
counterpart of the min-max problem, thus obtaining a new robust loss than can
be tractably minimized for many applications. Experimental results show that
this method has a performance edge for larger uncertainty sets. Both methods
run significantly faster than adversarial training, although their memory
requirements increase proportionally to the dimension of the input. Lastly, we
provide evidence that adding robustness can improve the natural accuracy of
neural networks for classification problems with tabular or vision data.

#### Code Availability Statement.

All the code to reproduce the results can be found here:
https://github.com/kimvc7/Robustness. An illustrative example on how to use
the code to train a network is provided here:
https://colab.research.google.com/github/kimvc7/Robustness/blob/main/demo.ipynb

#### Acknowledgements.

Xavier Boix has been supported by the Center for Brains, Minds and Machines
(funded by NSF STC award CCF- 1231216), the R01EY020517 grant from the
National Eye Institute (NIH) and Fujitsu Laboratories Ltd. (Contract No.
40008819).

## References

  * Anderson et al. (2020) Ross Anderson, Joey Huchette, Will Ma, Christian Tjandraatmadja, and Juan Pablo Vielma.  Strong mixed-integer programming formulations for trained neural networks.  _Mathematical Programming_ , pages 1–37, 2020.
  * Athalye et al. (2018) Anish Athalye, Nicholas Carlini, and David Wagner.  Obfuscated gradients give a false sense of security: Circumventing defenses to adversarial examples, 2018.
  * Balunovic and Vechev (2019) Mislav Balunovic and Martin Vechev.  Adversarial training and provable defenses: Bridging the gap.  In _International Conference on Learning Representations_ , 2019.
  * Ben-Tal et al. (2009) Aharon Ben-Tal, Laurent El Ghaoui, and Arkadi Nemirovski.  _Robust optimization_.  Princeton university press, 2009.
  * Ben-Tal et al. (2015) Aharon Ben-Tal, Dick Den Hertog, and Jean-Philippe Vial.  Deriving robust counterparts of nonlinear uncertain inequalities.  _Mathematical programming_ , 149(1):265–299, 2015\.
  * Bertsimas and den Hertog (2022) Dimitris Bertsimas and Dick den Hertog.  _Robust and Adaptive Optimization_.  Dynamic Ideas LLC, 2022.
  * Bertsimas et al. (2011) Dimitris Bertsimas, David B Brown, and Constantine Caramanis.  Theory and applications of robust optimization.  _SIAM Review_ , 53(3):464–501, 2011.
  * Bertsimas et al. (2019) Dimitris Bertsimas, Jack Dunn, Colin Pawlowski, and Ying Daisy Zhuo.  Robust classification.  _INFORMS Journal on Optimization_ , 1(1):2–34, 2019.
  * Bertsimas et al. (2020) Dimitris Bertsimas, Dick den Hertog, Jean Pauphilet, and Jianzhe Zhen.  Robust convex optimization: A new perspective that unifies and extends.  _Submitted to Mathematics of Operations Research_ , 2020.
  * Bunel et al. (2017) Rudy Bunel, Ilker Turkaslan, Philip HS Torr, Pushmeet Kohli, and M Pawan Kumar.  A unified view of piecewise linear neural network verification.  _arXiv preprint arXiv:1711.00455_ , 2017.
  * Chassein and Goerigk (2019) André Chassein and Marc Goerigk.  On the complexity of robust geometric programming with polyhedral uncertainty.  _Operations Research Letters_ , 47(1):21–24, 2019\.
  * Cohen et al. (2019) Jeremy Cohen, Elan Rosenfeld, and Zico Kolter.  Certified adversarial robustness via randomized smoothing.  In _International Conference on Machine Learning_ , pages 1310–1320. PMLR, 2019.
  * Dathathri et al. (2020) Sumanth Dathathri, Krishnamurthy Dvijotham, Alexey Kurakin, Aditi Raghunathan, Jonathan Uesato, Rudy Bunel, Shreya Shankar, Jacob Steinhardt, Ian Goodfellow, Percy Liang, and Pushmeet Kohli.  Enabling certification of verification-agnostic networks via memory-efficient semidefinite programming.  _arXiv preprint arXiv:2010.11645_ , 2020.
  * Deng (2012) Li Deng.  The mnist database of handwritten digit images for machine learning research.  _IEEE Signal Processing Magazine_ , 29(6):141–142, 2012.
  * Dua and Graff (2017) Dheeru Dua and Casey Graff.  UCI machine learning repository, 2017.  URL http://archive.ics.uci.edu/ml.
  * Dvijotham et al. (2018a) Krishnamurthy Dvijotham, Sven Gowal, Robert Stanforth, Relja Arandjelovic, Brendan O’Donoghue, Jonathan Uesato, and Pushmeet Kohli.  Training verified learners with learned verifiers.  _arXiv preprint arXiv:1805.10265_ , 2018a.
  * Dvijotham et al. (2018b) Krishnamurthy Dvijotham, Robert Stanforth, Sven Gowal, Timothy A Mann, and Pushmeet Kohli.  A dual approach to scalable verification of deep networks.  In _UAI_ , volume 1, page 3, 2018b.
  * Gehr et al. (2018) Timon Gehr, Matthew Mirman, Dana Drachsler-Cohen, Petar Tsankov, Swarat Chaudhuri, and Martin Vechev.  Ai2: Safety and robustness certification of neural networks with abstract interpretation.  In _2018 IEEE Symposium on Security and Privacy (SP)_ , pages 3–18. IEEE, 2018.
  * Goodfellow et al. (2015) Ian J. Goodfellow, Jonathon Shlens, and Christian Szegedy.  Explaining and harnessing adversarial examples, 2015.
  * Gowal et al. (2019) Sven Gowal, Krishnamurthy Dj Dvijotham, Robert Stanforth, Rudy Bunel, Chongli Qin, Jonathan Uesato, Relja Arandjelovic, Timothy Mann, and Pushmeet Kohli.  Scalable verified training for provably robust image classification.  In _Proceedings of the IEEE/CVF International Conference on Computer Vision_ , pages 4842–4851, 2019.
  * Hein and Andriushchenko (2017) Matthias Hein and Maksym Andriushchenko.  Formal guarantees on the robustness of a classifier against adversarial manipulation.  _CoRR_ , abs/1705.08475, 2017.  URL http://arxiv.org/abs/1705.08475.
  * Huang et al. (2016) Ruitong Huang, Bing Xu, Dale Schuurmans, and Csaba Szepesvari.  Learning with a strong adversary, 2016.
  * Ilyas et al. (2017) Andrew Ilyas, Ajil Jalal, Eirini Asteri, Constantinos Daskalakis, and Alexandros G. Dimakis.  The robust manifold defense: Adversarial training using generative models.  _CoRR_ , abs/1712.09196, 2017.  URL http://arxiv.org/abs/1712.09196.
  * Ilyas et al. (2019) Andrew Ilyas, Shibani Santurkar, Dimitris Tsipras, Logan Engstrom, Brandon Tran, and Aleksander Madry.  Adversarial examples are not bugs, they are features.  _arXiv preprint arXiv:1905.02175_ , 2019.
  * Kabilan et al. (2018) Vishaal Munusamy Kabilan, Brandon Morris, and Anh Nguyen.  Vectordefense: Vectorization as a defense to adversarial examples.  _CoRR_ , abs/1804.08529, 2018.  URL http://arxiv.org/abs/1804.08529.
  * Katz et al. (2017) Guy Katz, Clark Barrett, David L Dill, Kyle Julian, and Mykel J Kochenderfer.  Reluplex: An efficient smt solver for verifying deep neural networks.  In _International Conference on Computer Aided Verification_ , pages 97–117. Springer, 2017.
  * Krizhevsky (2009) Alex Krizhevsky.  Learning multiple layers of features from tiny images.  Technical report, 2009.
  * Kurakin et al. (2016) Alexey Kurakin, Ian Goodfellow, Samy Bengio, et al.  Adversarial examples in the physical world, 2016.
  * Lamb et al. (2018) Alex Lamb, Jonathan Binas, Anirudh Goyal, Dmitriy Serdyuk, Sandeep Subramanian, Ioannis Mitliagkas, and Yoshua Bengio.  Fortified networks: Improving the robustness of deep networks by modeling the manifold of hidden representations.  _arXiv preprint arXiv:1804.02485_ , 2018.
  * Lecuyer et al. (2019) Mathias Lecuyer, Vaggelis Atlidakis, Roxana Geambasu, Daniel Hsu, and Suman Jana.  Certified robustness to adversarial examples with differential privacy.  In _2019 IEEE Symposium on Security and Privacy (SP)_ , pages 656–672. IEEE, 2019.
  * Madry et al. (2019) Aleksander Madry, Aleksandar Makelov, Ludwig Schmidt, Dimitris Tsipras, and Adrian Vladu.  Towards deep learning models resistant to adversarial attacks, 2019.
  * Mirman et al. (2018) Matthew Mirman, Timon Gehr, and Martin Vechev.  Differentiable abstract interpretation for provably robust neural networks.  In _International Conference on Machine Learning_ , pages 3578–3586. PMLR, 2018.
  * Prakash et al. (2018) Aaditya Prakash, Nick Moran, Solomon Garber, Antonella DiLillo, and James A. Storer.  Deflecting adversarial attacks with pixel deflection.  _CoRR_ , abs/1801.08926, 2018.  URL http://arxiv.org/abs/1801.08926.
  * Raghunathan et al. (2018a) Aditi Raghunathan, Jacob Steinhardt, and Percy Liang.  Certified defenses against adversarial examples.  _arXiv preprint arXiv:1801.09344_ , 2018a.
  * Raghunathan et al. (2018b) Aditi Raghunathan, Jacob Steinhardt, and Percy Liang.  Semidefinite relaxations for certifying robustness to adversarial examples.  _arXiv preprint arXiv:1811.01057_ , 2018b.
  * Rauber et al. (2017) Jonas Rauber, Wieland Brendel, and Matthias Bethge.  Foolbox: A python toolbox to benchmark the robustness of machine learning models.  In _Reliable Machine Learning in the Wild Workshop, 34th International Conference on Machine Learning_ , 2017.  URL http://arxiv.org/abs/1707.04131.
  * Rauber et al. (2020) Jonas Rauber, Roland Zimmermann, Matthias Bethge, and Wieland Brendel.  Foolbox native: Fast adversarial attacks to benchmark the robustness of machine learning models in pytorch, tensorflow, and jax.  _Journal of Open Source Software_ , 5(53):2607, 2020.  doi: 10.21105/joss.02607.  URL https://doi.org/10.21105/joss.02607.
  * Rockafellar (1970) R. Tyrrell Rockafellar.  _Convex analysis_.  Princeton Mathematical Series. Princeton University Press, Princeton, N. J., 1970.
  * Roos et al. (2020) Ernst Roos, Dick den Hertog, Aharon Ben-Tal, F De Ruiter, and Jianzhe Zhen.  Tractable approximation of hard uncertain optimization problems.  _Available on Optimization-Online_ , 2020.
  * Ross and Doshi-Velez (2017) Andrew Slavin Ross and Finale Doshi-Velez.  Improving the adversarial robustness and interpretability of deep neural networks by regularizing their input gradients, 2017.
  * Singh et al. (2018) Gagandeep Singh, Timon Gehr, Matthew Mirman, Markus Püschel, and Martin T Vechev.  Fast and effective robustness certification.  _NeurIPS_ , 1(4):6, 2018.
  * Szegedy et al. (2013) Christian Szegedy, Wojciech Zaremba, Ilya Sutskever, Joan Bruna, Dumitru Erhan, Ian Goodfellow, and Rob Fergus.  Intriguing properties of neural networks.  _arXiv preprint arXiv:1312.6199_ , 2013.
  * Tjeng et al. (2019) Vincent Tjeng, Kai Xiao, and Russ Tedrake.  Evaluating robustness of neural networks with mixed integer programming, 2019.
  * Weng et al. (2018) Lily Weng, Huan Zhang, Hongge Chen, Zhao Song, Cho-Jui Hsieh, Luca Daniel, Duane Boning, and Inderjit Dhillon.  Towards fast computation of certified robustness for relu networks.  In _International Conference on Machine Learning_ , pages 5276–5285. PMLR, 2018.
  * Wong and Kolter (2018) Eric Wong and Zico Kolter.  Provable defenses against adversarial examples via the convex outer adversarial polytope.  In _International Conference on Machine Learning_ , pages 5286–5295. PMLR, 2018.
  * Wong et al. (2018) Eric Wong, Frank R Schmidt, Jan Hendrik Metzen, and J Zico Kolter.  Scaling provable adversarial defenses.  _arXiv preprint arXiv:1805.12514_ , 2018.
  * Xiao et al. (2017) Han Xiao, Kashif Rasul, and Roland Vollgraf.  Fashion-mnist: a novel image dataset for benchmarking machine learning algorithms, 2017.
  * Xie et al. (2017) Cihang Xie, Jianyu Wang, Zhishuai Zhang, Zhou Ren, and Alan Yuille.  Mitigating adversarial effects through randomization, 2017.
  * Yan et al. (2018) Ziang Yan, Yiwen Guo, and Changshui Zhang.  Deep defense: Training dnns with improved adversarial robustness.  _arXiv preprint arXiv:1803.00404_ , 2018.
  * Zhang et al. (2018) Huan Zhang, Tsui-Wei Weng, Pin-Yu Chen, Cho-Jui Hsieh, and Luca Daniel.  Efficient neural network robustness certification with general activation functions.  _arXiv preprint arXiv:1811.00866_ , 2018.
  * Zhang et al. (2019) Huan Zhang, Hongge Chen, Chaowei Xiao, Sven Gowal, Robert Stanforth, Bo Li, Duane Boning, and Cho-Jui Hsieh.  Towards stable and efficient training of verifiably robust neural networks.  _arXiv preprint arXiv:1906.06316_ , 2019.
  * Zhen et al. (2017) Jianzhe Zhen, FJCT de Ruiter, and D Den Hertog.  Robust optimization for models with uncertain soc and sdp constraints.  _Optimization Online_ , 2017.

## Appendix A: Proofs of Lemmas

In this Appendix, we prove lemmas 1, 2 and 3.

### A1: Proof of Lemma 1

Proof  Since $f$ is convex and closed, we have $f=(f^{*})^{*}$ (Rockafellar,
1970), and applying the definition of the convex conjugate function we obtain

$f({\boldsymbol{z}}(\boldsymbol{\delta}))=(f^{\star})^{\star}({\boldsymbol{z}}(\boldsymbol{\delta}))=\sup_{{\boldsymbol{u}}\in\text{dom}(f^{*})}{\boldsymbol{z}}(\boldsymbol{\delta})^{\top}{\boldsymbol{u}}-f^{\star}({\boldsymbol{u}}),$

which implies

$\displaystyle\sup_{\boldsymbol{\delta}\in\mathcal{U}}\>f(\boldsymbol{z}(\boldsymbol{\delta}))+g(\boldsymbol{z}(\boldsymbol{\delta}))$
$\displaystyle=\sup_{\boldsymbol{\delta}\in\mathcal{U}}\enspace\sup_{{\boldsymbol{u}}\in\text{dom}(f^{\star})}\>\boldsymbol{z}(\boldsymbol{\delta})^{T}{\boldsymbol{u}}-f^{\star}({\boldsymbol{u}})+g(\boldsymbol{z}(\boldsymbol{\delta}))$
$\displaystyle=\sup_{{\boldsymbol{u}}\in\text{dom}(f^{\star})}\enspace\sup_{\boldsymbol{\delta}\in\mathcal{U}}\>\boldsymbol{z}(\boldsymbol{\delta})^{T}{\boldsymbol{u}}-f^{\star}({\boldsymbol{u}})+g(\boldsymbol{z}(\boldsymbol{\delta})),$

as desired.

### A2: Proof of Lemma 2

Let
$\mathcal{Z}=\\{{\boldsymbol{z}}(\boldsymbol{\delta}):\boldsymbol{\delta}\in\mathcal{U}\\}$.
Defining the indicator function

$\gamma({\boldsymbol{z}}|\mathcal{Z})=\begin{cases}&0\quad\text{if
}{\boldsymbol{z}}\in\mathcal{Z},\\\ &\infty\quad\text{otherwise,}\end{cases}$

and applying the Fenchel duality theorem (Rockafellar, 1970), we obtain:

$\displaystyle\sup_{\boldsymbol{\delta}\in\mathcal{U}}g({\boldsymbol{z}}(\boldsymbol{\delta}))=\sup_{{\boldsymbol{z}}\in\mathcal{Z}}g({\boldsymbol{z}})=\sup_{{\boldsymbol{z}}\in\text{dom}(g)\cap\text{dom}(\gamma)}g({\boldsymbol{z}})-\gamma({\boldsymbol{z}}|\mathcal{Z})=\inf_{{\boldsymbol{v}}\in\text{dom}(g_{\star})}\gamma^{\star}({\boldsymbol{v}}|\mathcal{Z})-g_{\star}({\boldsymbol{v}})$
(28)

Finally, since
$\gamma^{\star}({\boldsymbol{v}}|\mathcal{Z})=\sup_{{\boldsymbol{z}}\in\mathcal{Z}}{\boldsymbol{z}}^{\top}{\boldsymbol{v}}$,
we conclude

$\sup_{\boldsymbol{\delta}\in\mathcal{U}}g({\boldsymbol{z}}(\boldsymbol{\delta}))=\inf_{{\boldsymbol{v}}\in\text{dom}(g_{\star})}\sup_{{\boldsymbol{z}}\in\mathcal{Z}}{\boldsymbol{z}}^{\top}{\boldsymbol{v}}-g_{\star}({\boldsymbol{v}})=\inf_{{\boldsymbol{v}}\in\text{dom}(g_{\star})}\sup_{\boldsymbol{\delta}\in\mathcal{U}}{\boldsymbol{z}}(\boldsymbol{\delta})^{\top}{\boldsymbol{v}}-g_{\star}({\boldsymbol{v}}).$

### A3: Proof of Lemma 3

We first prove part a). By definition, we have

$f^{\star}({\boldsymbol{z}})=\sup_{{\boldsymbol{x}}}{\boldsymbol{z}}^{\top}{\boldsymbol{x}}-{\boldsymbol{p}}^{\top}[{\boldsymbol{x}}]^{+}.$

Notice that if the $i^{th}$ component of ${\boldsymbol{z}}$ is negative for
any $i$, then $f^{\star}({\boldsymbol{z}})=\infty$ because ${\boldsymbol{x}}$
can be the vector with an arbitrarily large negative value in the $i^{th}$
coordinate and $0$ everywhere else. Similarly, if the $i^{th}$ component of
${\boldsymbol{z}}$ is larger than the $i^{th}$ coordinate of
${\boldsymbol{p}}$ for any $i$, then again
$f^{\star}({\boldsymbol{z}})=\infty$ because ${\boldsymbol{x}}$ can be the
vector with an arbitrarily large positive value in the $i^{th}$ coordinate and
$0$ everywhere else. Moreover, if
$\boldsymbol{0}\leq{\boldsymbol{z}}\leq{\boldsymbol{p}}$, then

$\sup_{{\boldsymbol{x}}}{\boldsymbol{z}}^{\top}{\boldsymbol{x}}-{\boldsymbol{p}}^{\top}[{\boldsymbol{x}}]^{+}\leq\sup_{{\boldsymbol{x}}}{\boldsymbol{z}}^{\top}{\boldsymbol{x}}-{\boldsymbol{z}}^{\top}[{\boldsymbol{x}}]^{+}=\sup_{{\boldsymbol{x}}}{\boldsymbol{z}}^{\top}({\boldsymbol{x}}-[{\boldsymbol{x}}]^{+})\leq
0.$

Since ${\boldsymbol{x}}=\boldsymbol{0}$ achieves an objective value of $0$, we
conclude that $\boldsymbol{0}\leq{\boldsymbol{z}}\leq{\boldsymbol{p}}$ implies
$f^{\star}({\boldsymbol{z}})=0$ as desired.

Next, we proceed to prove part b). By definition of the concave conjugate we
have

$g_{\star}({\boldsymbol{z}})=\inf_{{\boldsymbol{x}}}{\boldsymbol{z}}^{\top}{\boldsymbol{x}}-({\boldsymbol{x}}^{\top}{\boldsymbol{u}}-{\boldsymbol{q}}^{\top}[{\boldsymbol{x}}]^{+})=\inf_{{\boldsymbol{x}}}({\boldsymbol{z}}-{\boldsymbol{u}})^{\top}{\boldsymbol{x}}+{\boldsymbol{q}}^{\top}[{\boldsymbol{x}}]^{+}.$

If the $i^{th}$ component of ${\boldsymbol{z}}$ is larger than the $i^{th}$
component of ${\boldsymbol{u}}$ for any $i$, then
$g_{\star}({\boldsymbol{z}})=\infty$ because ${\boldsymbol{x}}$ can be the
vector with an arbitrarily large negative value in the $i^{th}$ coordinate and
$0$ everywhere else. Similarly, if the $i^{th}$ component of
${\boldsymbol{z}}$ is smaller than the $i^{th}$ coordinate of
$({\boldsymbol{u}}-{\boldsymbol{q}})$ for any $i$, then again
$g_{\star}({\boldsymbol{z}})=\infty$ because ${\boldsymbol{x}}$ can be the
vector with an arbitrarily large positive value in the $i^{th}$ coordinate and
$0$ everywhere else. In addition, if
$\boldsymbol{{\boldsymbol{u}}-{\boldsymbol{q}}}\leq{\boldsymbol{z}}\leq{\boldsymbol{u}}$,
then

$\inf_{{\boldsymbol{x}}}({\boldsymbol{z}}-{\boldsymbol{u}})^{\top}{\boldsymbol{x}}+{\boldsymbol{q}}^{\top}[{\boldsymbol{x}}]^{+}\geq\inf_{{\boldsymbol{x}}}({\boldsymbol{z}}-{\boldsymbol{u}})^{\top}{\boldsymbol{x}}+({\boldsymbol{u}}-{\boldsymbol{z}})^{\top}[{\boldsymbol{x}}]^{+}=\inf_{{\boldsymbol{x}}}({\boldsymbol{u}}-{\boldsymbol{z}})^{\top}([{\boldsymbol{x}}]^{+}-{\boldsymbol{x}})\geq
0.$

Since ${\boldsymbol{x}}=\boldsymbol{0}$ achieves an objective value of $0$, we
conclude that
$\boldsymbol{{\boldsymbol{u}}-{\boldsymbol{q}}}\leq{\boldsymbol{z}}\leq{\boldsymbol{u}}$
implies $g_{\star}({\boldsymbol{z}})=0$ as desired.

## Appendix B: Generalized Results

We now state and proof the generalization of Theorem 4, Corollary 5 and
Theorem 6 for the case in which the neural network has more than 2 layers.

######  Theorem 7 (Generalization of Theorem 4)

For all $2\leq l\leq L$, it holds

$\displaystyle\sup_{{\boldsymbol{\delta}}\in\mathcal{U}}\>$
$\displaystyle{\boldsymbol{c}}_{k}^{\top}{\boldsymbol{z}}^{\ell}(\theta,{\boldsymbol{x}}+{\boldsymbol{\delta}})=$
(29)
$\displaystyle\begin{split}&\sup_{{{\boldsymbol{s}}}_{L}}\inf_{{{\boldsymbol{t}}}_{L}}\dots\sup_{{{\boldsymbol{s}}}_{l}}\inf_{{{\boldsymbol{t}}}_{l}}\sup_{\boldsymbol{\delta}\in\mathcal{U}}\>({\boldsymbol{p}}_{l}-{\boldsymbol{q}}_{l})^{\top}{\boldsymbol{z}}^{l-1}(\theta,{\boldsymbol{x}}+{\boldsymbol{\delta}})+\sum_{\ell=l}^{L-1}({\boldsymbol{p}}_{\ell+1}-{\boldsymbol{q}}_{\ell+1})^{\top}{\boldsymbol{b}}^{\ell}+{\boldsymbol{c}}_{k}^{\top}{\boldsymbol{b}}^{L}\\\
&\hskip
99.58464pt\text{s.t.}\quad{\boldsymbol{p}}_{L}=[({\boldsymbol{W}}^{L})^{\top}{\boldsymbol{c}}_{k}]^{+}\odot{\boldsymbol{s}}_{L}\\\
&\hskip
125.19194pt{\boldsymbol{q}}_{L}=[-({\boldsymbol{W}}^{L})^{\top}{\boldsymbol{c}}_{k}]^{+}\odot{\boldsymbol{t}}_{L}\\\
&\hskip
125.19194pt{\boldsymbol{p}}_{\ell}=(([{\boldsymbol{W}}^{\ell}]^{+})^{\top}{\boldsymbol{p}}_{\ell+1}+([-{\boldsymbol{W}}^{\ell}]^{+})^{\top}{\boldsymbol{q}}_{\ell+1})\odot{\boldsymbol{s}}_{\ell}\quad\forall\>\ell=l,\dots,L-1\\\
&\hskip
125.19194pt{\boldsymbol{q}}_{\ell}=(([-{\boldsymbol{W}}^{\ell}]^{+})^{\top}{\boldsymbol{p}}_{\ell+1}+([{\boldsymbol{W}}^{\ell}]^{+})^{\top}{\boldsymbol{q}}_{\ell+1})\odot{\boldsymbol{t}}_{\ell}\quad\forall\>\ell=l,\dots,L-1\\\
&\hskip 125.19194pt0\leq{\boldsymbol{s}}_{\ell},{\boldsymbol{t}}_{\ell}\leq
1\quad\forall\>\ell=l,\dots,L.\end{split}$ (30)

Proof  We will proceed by backward induction on the layer number $l$.
Case $l=L$:
The proof is equivalent to the case $L=2$ already proved in Section 5.

Case $l-1$:
Suppose the theorem holds for some fixed $l$ with $l>2$. We have

$\displaystyle({\boldsymbol{p}}_{l}-{\boldsymbol{q}}_{l})^{\top}{\boldsymbol{z}}^{l-1}(\theta,{\boldsymbol{x}}+{\boldsymbol{\delta}})$
$\displaystyle=$
$\displaystyle({\boldsymbol{p}}_{l}-{\boldsymbol{q}}_{l})^{\top}({\boldsymbol{W}}^{l-1}[{\boldsymbol{z}}^{l-2}(\theta,{\boldsymbol{x}}+{\boldsymbol{\delta}})]^{+}+{\boldsymbol{b}}^{l-1})$
$\displaystyle=$ $\displaystyle
f_{+}({\boldsymbol{z}}^{l-2}(\theta,{\boldsymbol{x}}+{\boldsymbol{\delta}}))-f_{-}({\boldsymbol{z}}^{l-2}(\theta,{\boldsymbol{x}}+{\boldsymbol{\delta}}))+({\boldsymbol{p}}_{l}-{\boldsymbol{q}}_{l})^{\top}{\boldsymbol{b}}^{l-1},$

where

$\displaystyle f_{+}({\boldsymbol{x}})$
$\displaystyle=({\boldsymbol{p}}_{l}^{\top}[{\boldsymbol{W}}^{l-1}]^{+}+{\boldsymbol{q}}_{l}^{\top}[-{\boldsymbol{W}}^{l-1}]^{+})[{\boldsymbol{x}}]^{+},\quad\text{and}$
$\displaystyle f_{-}({\boldsymbol{x}})$
$\displaystyle=({\boldsymbol{p}}_{l}^{\top}[-{\boldsymbol{W}}^{l-1}]^{+}+{\boldsymbol{q}}_{l}^{\top}[{\boldsymbol{W}}^{l-1}]^{+})[{\boldsymbol{x}}]^{+}.$

By Lemma 1 we then obtain

$\displaystyle\sup_{{\boldsymbol{\delta}}\in\mathcal{U}}\>({\boldsymbol{p}}_{l}-{\boldsymbol{q}}_{l})^{\top}{\boldsymbol{z}}^{l-1}(\theta,{\boldsymbol{x}}+{\boldsymbol{\delta}})$
$\displaystyle=$
$\displaystyle\sup_{{\boldsymbol{\delta}}\in\mathcal{U}}\>f_{+}({\boldsymbol{z}}^{l-2}(\theta,{\boldsymbol{x}}+{\boldsymbol{\delta}}))-f_{-}({\boldsymbol{z}}^{l-2}(\theta,{\boldsymbol{x}}+{\boldsymbol{\delta}}))+({\boldsymbol{p}}_{l}-{\boldsymbol{q}}_{l})^{\top}{\boldsymbol{b}}^{l-1}$
$\displaystyle=$
$\displaystyle\sup_{{\boldsymbol{u}}_{l-1}\in\text{dom}(f_{+}^{\star})}\sup_{{\boldsymbol{\delta}}\in\mathcal{U}}{\boldsymbol{u}}_{l-1}^{\top}{\boldsymbol{z}}^{l-2}(\theta,{\boldsymbol{x}}+{\boldsymbol{\delta}})-f_{-}({\boldsymbol{z}}^{l-2}(\theta,{\boldsymbol{x}}+{\boldsymbol{\delta}}))+({\boldsymbol{p}}_{l}-{\boldsymbol{q}}_{l})^{\top}{\boldsymbol{b}}^{l-1}.$

Defining the concave function
$g({\boldsymbol{x}})={\boldsymbol{u}}_{l-1}^{\top}{\boldsymbol{x}}-f_{-}({\boldsymbol{x}})$,
and applying Lemma 2 we obtain

$\displaystyle\sup_{{\boldsymbol{\delta}}\in\mathcal{U}}\>({\boldsymbol{p}}_{l}-{\boldsymbol{q}}_{l})^{\top}{\boldsymbol{z}}^{l-1}(\theta,{\boldsymbol{x}}+{\boldsymbol{\delta}})$
(31) $\displaystyle=$
$\displaystyle\sup_{{\boldsymbol{u}}_{l-1}\in\text{dom}(f_{+}^{\star})}\inf_{{\boldsymbol{v}}_{l-1}\in\text{dom}(g_{\star})}\sup_{{\boldsymbol{\delta}}\in\mathcal{U}}{\boldsymbol{v}}_{l-1}^{\top}{\boldsymbol{z}}^{l-2}(\theta,{\boldsymbol{x}}+{\boldsymbol{\delta}})+({\boldsymbol{p}}_{l}-{\boldsymbol{q}}_{l})^{\top}{\boldsymbol{b}}^{l-1}.$
(32)

Lastly, by Lemma 3 we can substitute

$\displaystyle{\boldsymbol{u}}_{l-1}$
$\displaystyle=([({\boldsymbol{W}}^{l-1}]^{+})^{\top}{\boldsymbol{p}}_{l}+([-{\boldsymbol{W}}^{l-1}]^{+})^{\top}{\boldsymbol{q}}_{l})\odot{\boldsymbol{s}}_{l-1}$
$\displaystyle={\boldsymbol{p}}_{l-1},\quad\text{and}$
$\displaystyle{\boldsymbol{v}}_{l-1}$
$\displaystyle=(([{\boldsymbol{W}}^{l-1}]^{+})^{\top}{\boldsymbol{p}}_{l}+([-{\boldsymbol{W}}^{l-1}]^{+})^{\top}{\boldsymbol{q}}_{l})\odot{\boldsymbol{s}}_{l-1}-({\boldsymbol{p}}_{l}[-{\boldsymbol{W}}^{l-1}]^{+}+{\boldsymbol{q}}_{l}[{\boldsymbol{W}}^{l-1}]^{+})\odot{\boldsymbol{t}}_{l-1}$
$\displaystyle={\boldsymbol{p}}_{l-1}-{\boldsymbol{q}}_{l-1},$

which together with the induction hypothesis imply that Problem (30) is
equivalent to

$\displaystyle\sup_{{\boldsymbol{\delta}}\in\mathcal{U}}\>$
$\displaystyle{\boldsymbol{c}}_{k}^{\top}{\boldsymbol{z}}^{\ell-1}(\theta,{\boldsymbol{x}}+{\boldsymbol{\delta}})=$
$\displaystyle\sup_{{{\boldsymbol{s}}}_{L}}\inf_{{{\boldsymbol{t}}}_{L}}\dots\sup_{{{\boldsymbol{s}}}_{l-1}}\inf_{{{\boldsymbol{t}}}_{l-1}}\sup_{\boldsymbol{\delta}\in\mathcal{U}}\>({\boldsymbol{p}}_{l-1}-{\boldsymbol{q}}_{l-1})^{\top}{\boldsymbol{z}}^{l-2}(\theta,{\boldsymbol{x}}+{\boldsymbol{\delta}})+\sum_{\ell=l-1}^{L-1}({\boldsymbol{p}}_{\ell+1}-{\boldsymbol{q}}_{\ell+1})^{\top}{\boldsymbol{b}}^{\ell}+{\boldsymbol{c}}_{k}^{\top}{\boldsymbol{b}}^{L}$
$\displaystyle\hskip
99.58464pt\text{s.t.}\quad{\boldsymbol{p}}_{L}=[({\boldsymbol{W}}^{L})^{\top}{\boldsymbol{c}}_{k}]^{+}\odot{\boldsymbol{s}}_{L}$
$\displaystyle\hskip
125.19194pt{\boldsymbol{q}}_{L}=[-({\boldsymbol{W}}^{L})^{\top}{\boldsymbol{c}}_{k}]^{+}\odot{\boldsymbol{t}}_{L}$
$\displaystyle\hskip
125.19194pt{\boldsymbol{p}}_{\ell}=(([{\boldsymbol{W}}^{\ell}]^{+})^{\top}{\boldsymbol{p}}_{\ell+1}+([-{\boldsymbol{W}}^{\ell}]^{+})^{\top}{\boldsymbol{q}}_{\ell+1})\odot{\boldsymbol{s}}_{\ell}\quad\forall\>l-1\leq\ell\leq
L-1$ $\displaystyle\hskip
125.19194pt{\boldsymbol{q}}_{\ell}=(([-{\boldsymbol{W}}^{\ell}]^{+})^{\top}{\boldsymbol{p}}_{\ell+1}+([{\boldsymbol{W}}^{\ell}]^{+})^{\top}{\boldsymbol{q}}_{\ell+1})\odot{\boldsymbol{t}}_{\ell}\quad\forall\>l-1\leq\ell\leq
L-1$ $\displaystyle\hskip
125.19194pt0\leq{\boldsymbol{s}}_{\ell},{\boldsymbol{t}}_{\ell}\leq
1\quad\forall\>\ell=l-1,\dots,L,$

and therefore the theorem holds for $l-1$ as desired.

######  Corollary 8 (Generalization of Corollary 5)

If
$\mathcal{U}=\\{\boldsymbol{\delta}:\|\boldsymbol{\delta}\|_{p}\leq\rho\\}$,
then:

$\displaystyle\sup_{{\boldsymbol{\delta}}\in\mathcal{U}}\>{\boldsymbol{c}}_{k}^{\top}{\boldsymbol{z}}^{L}(\theta,{\boldsymbol{x}}+{\boldsymbol{\delta}})$
(33)
$\displaystyle\begin{split}=&\sup_{{{\boldsymbol{s}}}_{L}}\inf_{{{\boldsymbol{t}}}_{L}}\dots\sup_{{{\boldsymbol{s}}}_{2}}\inf_{{{\boldsymbol{t}}}_{2}}\rho\|({\boldsymbol{p}}_{2}-{\boldsymbol{q}}_{2})^{\top}{\boldsymbol{W}}^{1}\|_{q}+({\boldsymbol{p}}_{2}-{\boldsymbol{q}}_{2})^{\top}{\boldsymbol{W}}^{1}{\boldsymbol{x}}+\sum_{\ell=1}^{L-1}({\boldsymbol{p}}_{\ell+1}-{\boldsymbol{q}}_{\ell+1})^{\top}{\boldsymbol{b}}^{\ell}+{\boldsymbol{c}}_{k}^{\top}{\boldsymbol{b}}^{L}\\\
&\hskip
99.58464pt\text{s.t.}\quad{\boldsymbol{p}}_{L}=[({\boldsymbol{W}}^{L})^{\top}{\boldsymbol{c}}_{k}]^{+}\odot{\boldsymbol{s}}_{L}\\\
&\hskip
125.19194pt{\boldsymbol{q}}_{L}=[-({\boldsymbol{W}}^{L})^{\top}{\boldsymbol{c}}_{k}]^{+}\odot{\boldsymbol{t}}_{L}\\\
&\hskip
125.19194pt{\boldsymbol{p}}_{\ell}=(([{\boldsymbol{W}}^{\ell}]^{+})^{\top}{\boldsymbol{p}}_{\ell+1}+([-{\boldsymbol{W}}^{\ell}]^{+})^{\top}{\boldsymbol{q}}_{\ell+1})\odot{\boldsymbol{s}}_{\ell}\quad\forall\>\ell=2,\dots,L-1\\\
&\hskip
125.19194pt{\boldsymbol{q}}_{\ell}=(([-{\boldsymbol{W}}^{\ell}]^{+})^{\top}{\boldsymbol{p}}_{\ell+1}+([{\boldsymbol{W}}^{\ell}]^{+})^{\top}{\boldsymbol{q}}_{\ell+1})\odot{\boldsymbol{t}}_{\ell}\quad\forall\>\ell=2,\dots,L-1\\\
&\hskip 125.19194pt0\leq{\boldsymbol{s}}_{\ell},{\boldsymbol{t}}_{\ell}\leq
1\quad\forall\>\ell=2,\dots,L,\end{split}$ (34)

where $\|\cdot\|_{q}$ is the conjugate norm of $\|\cdot\|_{p}$.

Proof  The proof follows directly after applying Theorem 7 with $l=2$ and
using again Eq. (23).

###### Definition 9

We introduce the following definitions to simplify notation:

$\displaystyle{\boldsymbol{s}}$
$\displaystyle\coloneqq({\boldsymbol{s}}_{2},\ldots,{\boldsymbol{s}}_{L})$
$\displaystyle{\boldsymbol{t}}$
$\displaystyle\coloneqq({\boldsymbol{t}}_{2},\ldots,{\boldsymbol{t}}_{L})$
$\displaystyle{\boldsymbol{p}}_{L}({\boldsymbol{s}},{\boldsymbol{t}})$
$\displaystyle\coloneqq[({\boldsymbol{W}}^{2})^{\top}{\boldsymbol{c}}_{k}]^{+}\odot{\boldsymbol{s}}_{L}$
$\displaystyle{\boldsymbol{q}}_{L}({\boldsymbol{s}},{\boldsymbol{t}})$
$\displaystyle\coloneqq[-({\boldsymbol{W}}^{2})^{\top}{\boldsymbol{c}}_{k}]^{+}\odot{\boldsymbol{t}}_{L}$
$\displaystyle{\boldsymbol{p}}_{\ell}({\boldsymbol{s}},{\boldsymbol{t}})$
$\displaystyle\coloneqq\left(([{\boldsymbol{W}}^{\ell}]^{+})^{\top}{\boldsymbol{p}}_{\ell+1}({\boldsymbol{s}},{\boldsymbol{t}})+([-{\boldsymbol{W}}^{\ell}]^{+})^{\top}{\boldsymbol{q}}_{\ell+1}({\boldsymbol{s}},{\boldsymbol{t}})\right)\odot{\boldsymbol{s}}_{\ell}\quad\forall\>1\leq\ell<L$
$\displaystyle{\boldsymbol{q}}_{\ell}({\boldsymbol{s}},{\boldsymbol{t}})$
$\displaystyle\coloneqq\left(([-{\boldsymbol{W}}^{\ell}]^{+})^{\top}{\boldsymbol{p}}_{\ell+1}({\boldsymbol{s}},{\boldsymbol{t}})+([{\boldsymbol{W}}^{\ell}]^{+})^{\top}{\boldsymbol{q}}_{\ell+1}({\boldsymbol{s}},{\boldsymbol{t}})\right)\odot{\boldsymbol{t}}_{\ell}\quad\forall\>1\leq\ell<L$
$\displaystyle R_{\ell}({\boldsymbol{s}},{\boldsymbol{t}})$
$\displaystyle\coloneqq\sum_{\ell^{\prime}=\ell}^{L-1}\left({\boldsymbol{p}}_{\ell^{\prime}+1}({\boldsymbol{s}},{\boldsymbol{t}})-{\boldsymbol{q}}_{\ell^{\prime}+1}({\boldsymbol{s}},{\boldsymbol{t}})\right)^{\top}{\boldsymbol{b}}^{\ell^{\prime}}.$

######  Theorem 10 (Generalization of Theorem 6)

$\displaystyle\sup_{{\boldsymbol{\delta}}:\|\boldsymbol{\delta}\|_{1}\leq\rho}$
$\displaystyle\>{\boldsymbol{c}}_{k}^{\top}{\boldsymbol{z}}^{L}(\theta,{\boldsymbol{x}}+{\boldsymbol{\delta}})$
$\displaystyle\leq$ $\displaystyle\inf_{0\leq{{\boldsymbol{t}}}\leq
1}\>\max_{m\in[M]}\>\max\bigg{\\{}\
g^{L}_{k,m}(\theta,{\boldsymbol{x}},{\boldsymbol{t}},\rho),g^{L}_{k,m}(\theta,{\boldsymbol{x}},{\boldsymbol{t}},-\rho)\bigg{\\}},$
(35)

where the new network $g$ is defined by the equations

$\displaystyle
g^{1}_{m}({\boldsymbol{W}},{\boldsymbol{x}},{\boldsymbol{t}},a,r)$
$\displaystyle=r(a{\boldsymbol{W}}^{1}_{m}+{\boldsymbol{W}}^{1}{\boldsymbol{x}}+{\boldsymbol{b}}^{1})$
$\displaystyle g^{\ell}_{m}(\theta,{\boldsymbol{x}},{\boldsymbol{t}},a,r)$
$\displaystyle=[r{\boldsymbol{W}}^{\ell}]^{+}[g^{\ell-1}_{m}({\boldsymbol{W}},{\boldsymbol{x}},{\boldsymbol{t}},a,1)]^{+}+[-r{\boldsymbol{W}}^{\ell}]^{+}[g^{\ell-1}_{m}({\boldsymbol{W}},{\boldsymbol{x}},{\boldsymbol{t}},a,-1)]\odot{\boldsymbol{t}}_{\ell}+r{\boldsymbol{b}}^{\ell}$
$\displaystyle g^{L}_{k,m}(\theta,{\boldsymbol{x}},{\boldsymbol{t}},a)$
$\displaystyle=[{\boldsymbol{c}}_{k}^{\top}{\boldsymbol{W}}^{L}]^{+}[g^{L-1}_{m}(\theta,{\boldsymbol{x}},{\boldsymbol{t}},a,1)]^{+}+[-{\boldsymbol{c}}_{k}^{\top}{\boldsymbol{W}}^{L}]^{+}[g^{L-1}_{m}(\theta,{\boldsymbol{x}},{\boldsymbol{t}},a,-1)]\odot{\boldsymbol{t}}_{L}+{\boldsymbol{c}}_{k}^{\top}{\boldsymbol{b}}^{L},$

for all $1<\ell<L$, $1\leq k\leq K$, $a\in\\{\rho,-\rho\\}$, and
$r\in\\{-1,1\\}$.

The proof of this theorem relies on the following lemma.

###### Lemma 11

For all $2\leq\ell\leq L-1$ it holds

$\displaystyle\sup_{0\leq{{\boldsymbol{s}}}_{\ell}\leq
1}\>{\boldsymbol{p}}_{\ell}({\boldsymbol{s}},{\boldsymbol{t}})^{\top}g^{\ell-1}_{m}(\theta,{\boldsymbol{x}},{\boldsymbol{t}},a,1)+{\boldsymbol{q}}_{\ell}({\boldsymbol{s}},{\boldsymbol{t}})^{\top}g^{\ell-1}_{m}(\theta,{\boldsymbol{x}},{\boldsymbol{t}},a,-1)$
(36) $\displaystyle=$
$\displaystyle{\boldsymbol{p}}_{\ell+1}({\boldsymbol{s}},{\boldsymbol{t}})^{\top}g^{\ell}_{m}(\theta,{\boldsymbol{x}},{\boldsymbol{t}},a,1)+{\boldsymbol{q}}_{\ell+1}({\boldsymbol{s}},{\boldsymbol{t}})^{\top}g^{\ell}_{m}(\theta,{\boldsymbol{x}},{\boldsymbol{t}},a,-1)-({\boldsymbol{p}}_{\ell+1}({\boldsymbol{s}},{\boldsymbol{t}})-{\boldsymbol{q}}_{\ell+1}({\boldsymbol{s}},{\boldsymbol{t}}))^{\top}{\boldsymbol{b}}^{\ell}.$
(37)

Proof  Let $2\leq\ell\leq L-1$, we have

$\displaystyle\sup_{0\leq{{\boldsymbol{s}}}_{\ell}\leq
1}\>{\boldsymbol{p}}_{\ell}({\boldsymbol{s}},{\boldsymbol{t}})^{\top}g^{\ell-1}_{m}(\theta,{\boldsymbol{x}},{\boldsymbol{t}},a,1)+{\boldsymbol{q}}_{\ell}({\boldsymbol{s}},{\boldsymbol{t}})^{\top}g^{\ell-1}_{m}(\theta,{\boldsymbol{x}},{\boldsymbol{t}},a,-1)$
$\displaystyle=$ $\displaystyle\sup_{0\leq{{\boldsymbol{s}}}_{\ell}\leq
1}\>\left((([{\boldsymbol{W}}^{\ell}]^{+})^{\top}{\boldsymbol{p}}_{\ell+1}({\boldsymbol{s}},{\boldsymbol{t}})+([-{\boldsymbol{W}}^{\ell}]^{+})^{\top}{\boldsymbol{q}}_{\ell+1}({\boldsymbol{s}},{\boldsymbol{t}}))\odot{\boldsymbol{s}}_{\ell}\right)^{\top}g^{\ell-1}_{m}(\theta,{\boldsymbol{x}},{\boldsymbol{t}},a,1)+$
$\displaystyle\quad\quad\quad\left((([-{\boldsymbol{W}}^{\ell}]^{+})^{\top}{\boldsymbol{p}}_{\ell+1}({\boldsymbol{s}},{\boldsymbol{t}})+([{\boldsymbol{W}}^{\ell}]^{+})^{\top}{\boldsymbol{q}}_{\ell+1}({\boldsymbol{s}},{\boldsymbol{t}}))\odot{\boldsymbol{t}}_{\ell}\right)^{\top}g^{\ell-1}_{m}(\theta,{\boldsymbol{x}},{\boldsymbol{t}},a,-1)$
$\displaystyle=$ $\displaystyle\sup_{0\leq{{\boldsymbol{s}}}_{\ell}\leq
1}\>\left(([{\boldsymbol{W}}^{\ell}]^{+})^{\top}{\boldsymbol{p}}_{\ell+1}({\boldsymbol{s}},{\boldsymbol{t}})+([-{\boldsymbol{W}}^{\ell}]^{+})^{\top}{\boldsymbol{q}}_{\ell+1}({\boldsymbol{s}},{\boldsymbol{t}})\right)^{\top}\left(g^{\ell-1}_{m}(\theta,{\boldsymbol{x}},{\boldsymbol{t}},a,1)\odot{\boldsymbol{s}}_{\ell}\right)+$
$\displaystyle\quad\quad\quad\left(([-{\boldsymbol{W}}^{\ell}]^{+})^{\top}{\boldsymbol{p}}_{\ell+1}({\boldsymbol{s}},{\boldsymbol{t}})+([{\boldsymbol{W}}^{\ell}]^{+})^{\top}{\boldsymbol{q}}_{\ell+1}({\boldsymbol{s}},{\boldsymbol{t}})\right)^{\top}\left(g^{\ell-1}_{m}(\theta,{\boldsymbol{x}},{\boldsymbol{t}},a,-1)\odot{\boldsymbol{t}}_{\ell}\right)$
$\displaystyle=$
$\displaystyle\left(([{\boldsymbol{W}}^{\ell}]^{+})^{\top}{\boldsymbol{p}}_{\ell+1}({\boldsymbol{s}},{\boldsymbol{t}})+([-{\boldsymbol{W}}^{\ell}]^{+})^{\top}{\boldsymbol{q}}_{\ell+1}({\boldsymbol{s}},{\boldsymbol{t}})\right)^{\top}[g^{\ell-1}_{m}(\theta,{\boldsymbol{x}},{\boldsymbol{t}},a,1)]^{+}+$
$\displaystyle\left(([-{\boldsymbol{W}}^{\ell}]^{+})^{\top}{\boldsymbol{p}}_{\ell+1}({\boldsymbol{s}},{\boldsymbol{t}})+([{\boldsymbol{W}}^{\ell}]^{+})^{\top}{\boldsymbol{q}}_{\ell+1}({\boldsymbol{s}},{\boldsymbol{t}})\right)^{\top}\left(g^{\ell-1}_{m}(\theta,{\boldsymbol{x}},{\boldsymbol{t}},a,-1)\odot{\boldsymbol{t}}_{\ell}\right)$
$\displaystyle=$
$\displaystyle{\boldsymbol{p}}_{\ell+1}({\boldsymbol{s}},{\boldsymbol{t}})^{\top}g^{\ell}_{m}(\theta,{\boldsymbol{x}},{\boldsymbol{t}},a,1)+{\boldsymbol{q}}_{\ell+1}({\boldsymbol{s}},{\boldsymbol{t}})^{\top}g^{\ell}_{m}(\theta,{\boldsymbol{x}},{\boldsymbol{t}},a,-1)-({\boldsymbol{p}}_{\ell+1}({\boldsymbol{s}},{\boldsymbol{t}})-{\boldsymbol{q}}_{\ell+1}({\boldsymbol{s}},{\boldsymbol{t}}))^{\top}{\boldsymbol{b}}^{\ell},$

as desired.

Proof [Proof of theorem 10]
By Corollary 8 with $p=1$ we know

$\displaystyle\sup_{{\boldsymbol{\delta}}:\|{\boldsymbol{\delta}}\|_{1}\leq\rho}\>{\boldsymbol{c}}_{k}^{\top}({\boldsymbol{z}}^{L}(\theta,{\boldsymbol{x}}+{\boldsymbol{\delta}})-{\boldsymbol{b}}^{L})$
(38) $\displaystyle\leq$
$\displaystyle\sup_{{{\boldsymbol{s}}}_{L}}\inf_{{\boldsymbol{t}}_{L}},\dots\sup_{{\boldsymbol{s}}_{2}}\inf_{{{\boldsymbol{t}}}_{2}}\>\rho\|({\boldsymbol{p}}_{2}({\boldsymbol{s}},{\boldsymbol{t}})-{\boldsymbol{q}}_{2}({\boldsymbol{s}},{\boldsymbol{t}}))^{\top}{\boldsymbol{W}}^{1}\|_{\infty}+({\boldsymbol{p}}_{2}({\boldsymbol{s}},{\boldsymbol{t}})-{\boldsymbol{q}}_{2}({\boldsymbol{s}},{\boldsymbol{t}}))^{\top}{\boldsymbol{W}}^{1}{\boldsymbol{x}}+R_{1}({\boldsymbol{s}},{\boldsymbol{t}})$
(39)
$\displaystyle\begin{split}\leq&\inf_{{{\boldsymbol{t}}}_{L}\dots{\boldsymbol{t}}_{2}}\sup_{{{\boldsymbol{s}}}_{L}\dots{\boldsymbol{s}}_{2}}\rho\|({\boldsymbol{p}}_{2}({\boldsymbol{s}},{\boldsymbol{t}})-{\boldsymbol{q}}_{2}({\boldsymbol{s}},{\boldsymbol{t}}))^{\top}{\boldsymbol{W}}^{1}\|_{\infty}+({\boldsymbol{p}}_{2}({\boldsymbol{s}},{\boldsymbol{t}})-{\boldsymbol{q}}_{2}({\boldsymbol{s}},{\boldsymbol{t}}))^{\top}{\boldsymbol{W}}^{1}{\boldsymbol{x}}+R_{1}({\boldsymbol{s}},{\boldsymbol{t}}),\end{split}$
(40)

where the last inequality follows from the min-max inequality. Observe that
${\boldsymbol{p}}_{\ell^{\prime}}({\boldsymbol{s}},{\boldsymbol{t}})$ and
${\boldsymbol{q}}_{\ell^{\prime}}({\boldsymbol{s}},{\boldsymbol{t}})$ are
independent on ${\boldsymbol{s}}_{\ell}$ for all ${\ell}^{\prime}>\ell$, which
in turn implies that $R_{\ell^{\prime}}({\boldsymbol{s}},{\boldsymbol{t}})$
does not depend on ${\boldsymbol{s}}_{\ell}$ for all ${\ell}^{\prime}>\ell$.
We can then solve the optimization problem in Eq. (34) for fixed
${{\boldsymbol{t}}}$ as follows:

$\displaystyle\sup_{0\leq{\boldsymbol{s}}_{L},\dots,{{\boldsymbol{s}}}_{2}\leq
1}\>\rho\|({\boldsymbol{p}}_{2}({\boldsymbol{s}},{\boldsymbol{t}})-{\boldsymbol{q}}_{2}({\boldsymbol{s}},{\boldsymbol{t}}))^{\top}{\boldsymbol{W}}^{1}\|_{\infty}+({\boldsymbol{p}}_{2}({\boldsymbol{s}},{\boldsymbol{t}})-{\boldsymbol{q}}_{2}({\boldsymbol{s}},{\boldsymbol{t}}))^{\top}{\boldsymbol{W}}^{1}{\boldsymbol{x}}+R_{1}({\boldsymbol{s}},{\boldsymbol{t}})$
$\displaystyle=$
$\displaystyle\max_{m\in[M]}\>\max\bigg{\\{}\sup_{0\leq{\boldsymbol{s}}_{L},\dots,{{\boldsymbol{s}}}_{2}\leq
1}\>({\boldsymbol{p}}_{2}({\boldsymbol{s}},{\boldsymbol{t}})-{\boldsymbol{q}}_{2}({\boldsymbol{s}},{\boldsymbol{t}}))^{\top}({\boldsymbol{W}}^{1}({\boldsymbol{x}}+\rho{\boldsymbol{e}}_{m})+{\boldsymbol{b}}^{1})+R_{2}({\boldsymbol{s}},{\boldsymbol{t}}),$
$\displaystyle\hskip
59.75095pt\sup_{0\leq{\boldsymbol{s}}_{L},\dots,{{\boldsymbol{s}}}_{2}\leq
1}\>({\boldsymbol{p}}_{2}({\boldsymbol{s}},{\boldsymbol{t}})-{\boldsymbol{q}}_{2}({\boldsymbol{s}},{\boldsymbol{t}}))^{\top}({\boldsymbol{W}}^{1}({\boldsymbol{x}}-\rho{\boldsymbol{e}}_{m})+{\boldsymbol{b}}^{1})+R_{2}({\boldsymbol{s}},{\boldsymbol{t}})\bigg{\\}}$
$\displaystyle=$
$\displaystyle\max_{m\in[M]}\>\max\bigg{\\{}\sup_{0\leq{\boldsymbol{s}}_{L},\dots,{{\boldsymbol{s}}}_{2}\leq
1}\>{\boldsymbol{p}}_{2}({\boldsymbol{s}},{\boldsymbol{t}})^{\top}g^{1}_{m}(\theta,{\boldsymbol{x}},{\boldsymbol{t}},\rho,1)+{\boldsymbol{q}}_{2}({\boldsymbol{s}},{\boldsymbol{t}})^{\top}g^{1}_{m}(\theta,{\boldsymbol{x}},{\boldsymbol{t}},\rho,-1)+R_{2}({\boldsymbol{s}},{\boldsymbol{t}}),$
$\displaystyle\hskip
59.75095pt\sup_{0\leq{\boldsymbol{s}}_{L},\dots,{{\boldsymbol{s}}}_{2}\leq
1}\>{\boldsymbol{p}}_{2}({\boldsymbol{s}},{\boldsymbol{t}})^{\top}g^{1}_{m}(\theta,{\boldsymbol{x}},{\boldsymbol{t}},-\rho,1)+{\boldsymbol{q}}_{2}({\boldsymbol{s}},{\boldsymbol{t}})^{\top}g^{1}_{m}(\theta,{\boldsymbol{x}},{\boldsymbol{t}},-\rho,-1)+R_{2}({\boldsymbol{s}},{\boldsymbol{t}})\bigg{\\}}.$

By repeatedly applying Lemma 11 for each $\ell=2,\dots,L-1$ we obtain

$\displaystyle\sup_{0\leq{\boldsymbol{s}}_{L},\dots,{{\boldsymbol{s}}}_{2}\leq
1}\>\rho\|({\boldsymbol{p}}_{2}({\boldsymbol{s}},{\boldsymbol{t}})-{\boldsymbol{q}}_{2}({\boldsymbol{s}},{\boldsymbol{t}}))^{\top}{\boldsymbol{W}}^{1}\|_{\infty}+({\boldsymbol{p}}_{2}({\boldsymbol{s}},{\boldsymbol{t}})-{\boldsymbol{q}}_{2}({\boldsymbol{s}},{\boldsymbol{t}}))^{\top}{\boldsymbol{W}}^{1}{\boldsymbol{x}}+R_{1}({\boldsymbol{s}},{\boldsymbol{t}})$
$\displaystyle=$
$\displaystyle\max_{m\in[M]}\>\max\bigg{\\{}\sup_{0\leq{\boldsymbol{s}}_{L}\leq
1}\>{\boldsymbol{p}}_{L}({\boldsymbol{s}},{\boldsymbol{t}})^{\top}g^{L-1}_{m}(\theta,{\boldsymbol{x}},{\boldsymbol{t}},\rho,1)+{\boldsymbol{q}}_{L}({\boldsymbol{s}},{\boldsymbol{t}})^{\top}g^{L-1}_{m}(\theta,{\boldsymbol{x}},{\boldsymbol{t}},\rho,-1),$
$\displaystyle\hskip 59.75095pt\sup_{0\leq{\boldsymbol{s}}_{L}\leq
1}\>{\boldsymbol{p}}_{L}({\boldsymbol{s}},{\boldsymbol{t}})^{\top}g^{L-1}_{m}(\theta,{\boldsymbol{x}},{\boldsymbol{t}},-\rho,1)+{\boldsymbol{q}}_{L}({\boldsymbol{s}},{\boldsymbol{t}})^{\top}g^{L-1}_{m}(\theta,{\boldsymbol{x}},{\boldsymbol{t}},-\rho,-1)\bigg{\\}}$
$\displaystyle=$
$\displaystyle\max_{m\in[M]}\>\max\bigg{\\{}\sup_{0\leq{\boldsymbol{s}}_{L}\leq
1}[{\boldsymbol{c}}_{k}^{\top}{\boldsymbol{W}}^{L}]^{+}\\!\\!\left(g^{L-1}_{m}(\theta,{\boldsymbol{x}},{\boldsymbol{t}},\rho,1)\odot{\boldsymbol{s}}_{L}\right)\\!+\\![-{\boldsymbol{c}}_{k}^{\top}{\boldsymbol{W}}^{L}]^{+}\\!(g^{L-1}_{m}(\theta,{\boldsymbol{x}},{\boldsymbol{t}},\rho,-1)\odot{\boldsymbol{t}}_{L}),$
$\displaystyle\hskip 59.75095pt\sup_{0\leq{\boldsymbol{s}}_{L}\leq
1}[{\boldsymbol{c}}_{k}^{\top}{\boldsymbol{W}}^{L}]^{+}\\!\\!\left(g^{L-1}_{m}(\theta,{\boldsymbol{x}},{\boldsymbol{t}},-\rho,1)\odot{\boldsymbol{s}}_{L}\right)\\!+\\![-{\boldsymbol{c}}_{k}^{\top}{\boldsymbol{W}}^{L}]^{+}(g^{L-1}_{m}(\theta,{\boldsymbol{x}},{\boldsymbol{t}},-\rho,-1)\odot{\boldsymbol{t}}_{L})\\!\bigg{\\}}$
$\displaystyle=$
$\displaystyle\max_{m\in[M]}\>\max\bigg{\\{}g^{L}_{k,m}(\theta,{\boldsymbol{x}},{\boldsymbol{t}},\rho),g^{L}_{k,m}(\theta,{\boldsymbol{x}},{\boldsymbol{t}},-\rho)\bigg{\\}}-{\boldsymbol{c}}_{k}^{\top}{\boldsymbol{b}}^{L},$

as desired.

## Appendix C: Other Adversarial Attacks for Vision Data Sets

In this section we further analyze adversarial accuracies with respect to the
vision data sets by evaluating other types of attacks. More specifically, we
show the adversarial accuracies of each method with respect to PGD attacks
bouned in $L_{1}$ and $L_{\infty}$ norms. Accross the three data sets (Fashion
MNIST, MNIST and CIFAR), we observe similar trends as those obtained with
$L_{2}$ norm bounded attacks in Section 6.3.

$\boldsymbol{\rho}=$ | 0.00 | 0.01 | 0.06 | 0.28 | 2.80 | 8.40 | 14.00 | 28.00
---|---|---|---|---|---|---|---|---
RUB | 90.27 | 90.20 | 90.16 | 90.27 | 88.36 | 86.25 | 84.18 | 76.25
ARUB-$L_{1}$ | 90.51 | 90.51 | 90.43 | 90.35 | 89.02 | 86.52 | 84.34 | 79.41
ARUB-$L_{\infty}$ | 89.96 | 89.96 | 89.92 | 89.84 | 88.44 | 84.69 | 81.41 | 74.30
Baseline-$L_{\infty}$ | 89.49 | 89.49 | 90.31 | 89.80 | 87.93 | 84.26 | 80.23 | 70.94
PGD-$L_{\infty}$ | 89.92 | 89.84 | 90.78 | 90.55 | 88.52 | 85.35 | 83.95 | 75.20
Nominal | 88.59 | 88.59 | 88.55 | 88.52 | 87.50 | 84.02 | 80.70 | 71.33
Table 9: Adversarial Accuracy for Fashion MNIST with PGD-$L_{1}$ attacks. $\boldsymbol{\rho}=$ | 0.00 | 0.01 | 0.06 | 0.28 | 2.80 | 8.40 | 14.00 | 28.00
---|---|---|---|---|---|---|---|---
RUB | 98.01 | 97.97 | 97.93 | 97.89 | 97.30 | 97.19 | 96.25 | 93.32
ARUB-$L_{1}$ | 98.52 | 98.52 | 98.48 | 98.48 | 98.05 | 97.50 | 96.52 | 94.92
ARUB-$L_{\infty}$ | 98.40 | 98.40 | 98.40 | 98.24 | 98.32 | 97.85 | 97.11 | 94.53
Baseline-$L_{\infty}$ | 98.05 | 98.05 | 98.05 | 98.05 | 97.58 | 96.52 | 95.55 | 92.30
PGD-$L_{\infty}$ | 98.48 | 98.48 | 98.48 | 98.44 | 98.28 | 97.73 | 97.30 | 95.00
Nominal | 97.73 | 97.73 | 97.73 | 97.54 | 97.30 | 95.94 | 94.49 | 90.94
Table 10: Adversarial Accuracy for MNIST with PGD-$L_{1}$ attacks. $\boldsymbol{\rho}=$ | 0.00 | 0.01 | 0.06 | 0.55 | 5.54 | 16.63 | 27.71 | 55.43
---|---|---|---|---|---|---|---|---
RUB | 53.59 | 53.59 | 53.59 | 50.94 | 51.33 | 48.28 | 45.94 | 41.33
ARUB-$L_{1}$ | 53.40 | 53.32 | 53.09 | 51.95 | 50.90 | 48.87 | 46.60 | 41.37
ARUB-$L_{\infty}$ | 53.67 | 53.63 | 53.52 | 52.97 | 45.51 | 42.19 | 40.62 | 37.38
Baseline-$L_{\infty}$ | 53.32 | 53.28 | 53.24 | 51.17 | 43.01 | 41.29 | 37.77 | 28.83
PGD-$L_{\infty}$ | 54.88 | 54.88 | 54.80 | 53.71 | 48.91 | 47.15 | 45.70 | 41.60
Nominal | 46.88 | 46.88 | 46.88 | 46.52 | 43.91 | 40.51 | 37.46 | 29.73
Table 11: Adversarial Accuracy for CIFAR with PGD-$L_{1}$ attacks. $\boldsymbol{\rho}=$ | 0.000 | 0.001 | 0.003 | 0.010 | 0.100 | 0.300 | 0.500 | 1.000
---|---|---|---|---|---|---|---|---
RUB | 89.38 | 89.18 | 88.71 | 87.15 | 66.48 | 56.17 | 56.13 | 55.78
ARUB-$L_{1}$ | 90.00 | 89.77 | 89.22 | 87.27 | 67.07 | 31.05 | 11.02 | 9.84
ARUB-$L_{\infty}$ | 89.77 | 89.22 | 88.87 | 86.76 | 80.43 | 67.93 | 54.49 | 12.15
Baseline-$L_{\infty}$ | 90.00 | 88.79 | 87.93 | 85.59 | 56.80 | 30.20 | 29.06 | 27.19
PGD-$L_{\infty}$ | 89.02 | 89.18 | 88.24 | 87.77 | 80.82 | 73.09 | 65.74 | 31.21
Nominal | 88.16 | 87.97 | 87.30 | 86.05 | 58.32 | 21.80 | 18.79 | 16.37
Table 12: Adversarial Accuracy for Fashion MNIST with PGD-$L_{\infty}$ attacks. $\boldsymbol{\rho}=$ | 0.000 | 0.001 | 0.003 | 0.010 | 0.100 | 0.300 | 0.500 | 1.000
---|---|---|---|---|---|---|---|---
RUB | 98.52 | 98.36 | 98.32 | 98.36 | 86.09 | 68.52 | 67.15 | 61.56
ARUB-$L_{1}$ | 99.14 | 99.02 | 98.87 | 98.44 | 87.97 | 61.88 | 59.61 | 54.77
ARUB-$L_{\infty}$ | 98.91 | 98.83 | 98.71 | 98.79 | 97.42 | 91.52 | 80.98 | 49.69
Baseline-$L_{\infty}$ | 98.59 | 98.09 | 98.63 | 97.66 | 84.69 | 65.08 | 63.24 | 59.49
PGD-$L_{\infty}$ | 98.71 | 98.87 | 98.83 | 98.75 | 97.66 | 93.83 | 86.80 | 64.57
Nominal | 98.05 | 97.97 | 97.85 | 97.30 | 81.91 | 23.28 | 18.98 | 17.38
Table 13: Adversarial Accuracy for MNIST with PGD-$L_{\infty}$ attacks. $\boldsymbol{\rho}=$ | 0.000 | 0.001 | 0.003 | 0.010 | 0.100 | 0.300 | 0.500 | 1.000
---|---|---|---|---|---|---|---|---
RUB | 50.86 | 46.02 | 48.05 | 45.94 | 21.64 | 10.62 | 10.62 | 10.62
ARUB-$L_{1}$ | 53.28 | 51.48 | 50.55 | 47.77 | 27.19 | 10.66 | 10.12 | 10.12
ARUB-$L_{\infty}$ | 55.43 | 51.21 | 45.66 | 42.70 | 29.61 | 10.43 | 11.17 | 11.05
Baseline-$L_{\infty}$ | 51.48 | 47.54 | 41.56 | 39.22 | 9.92 | 9.92 | 9.92 | 9.92
PGD-$L_{\infty}$ | 54.49 | 51.60 | 48.05 | 46.33 | 25.31 | 16.09 | 12.93 | 12.89
Nominal | 49.26 | 46.17 | 43.44 | 39.18 | 9.92 | 9.92 | 9.92 | 9.92
Table 14: Adversarial Accuracy for CIFAR with PGD-$L_{\infty}$ attacks.

## Appendix D: Convolutional Neural Networks

While in this paper we only consider feed forward neural networks, it is
possible to extend our results to convolutional neural networks that use ReLU
and MaxPool activation functions. In fact, Lemma 3 can be modified as follows:

###### Lemma 12

Let $x\in\mathbb{R}^{A\times B}$. Define $MP(x)\in\mathbb{R}^{C\times D}$ as
the MaxPool function whose $(c,d)$ coordinate corresponds to $\max_{i\in
I_{cd}}\\{x_{i}\\}$ for fixed sets of indices $I_{cd}$, and denote by
$\circledast$ the convolution operation. If
${\boldsymbol{u}}\in\mathbb{R}^{A\times B}$,
${\boldsymbol{p}},{\boldsymbol{q}}\in\mathbb{R}^{C\times D}$ have all
nonnegative coordinates, then the functions
$f({\boldsymbol{x}})={\boldsymbol{p}}\circledast MP[{\boldsymbol{x}}]^{+}$ and
$g({\boldsymbol{x}})={\boldsymbol{x}}\circledast{\boldsymbol{u}}-{\boldsymbol{q}}\circledast
MP[{\boldsymbol{x}}]^{+}$ satisfy

$\displaystyle a)\>f^{\star}({\boldsymbol{z}})=\begin{cases}0&\quad\text{if
}0\leq\sum_{i\in I_{cd}}z_{i}\leq p_{cd},\forall\>c\in[C],d\in[D],\\\
\infty&\quad\text{otherwise,}\end{cases}\quad\text{and}$ $\displaystyle
b)\>g_{\star}({\boldsymbol{z}})=\begin{cases}0&\quad\text{if
}u_{cd}-q_{cd}\leq\sum_{i\in I_{cd}}z_{i}\leq
u_{cd},\forall\>c\in[C],d\in[D],\\\ -\infty&\quad\text{otherwise.}\end{cases}$

The lemma above allows to obtain an upper bound for the robust loss of
convolutional networks with a very similar proof to that of Theorem 10.
However, convolutional networks are notoriously more memory consuming and
therefore computation of the robust upper bound requires more resources. We
have then left this computation for future work, but we here report results
for the other robust training methods using these more complex neural
networks.

We evaluate a convolutional neural network (denoted as _CNN_) that has been
commonly used in previous works of adversarial robustness Madry et al. (2019).
It has two convolutional layers alternated with pooling operations, and two
dense layers. We compare adversarial accuracy accross four different methods:
ARUB-$L_{\infty}$, ARUB-$L_{1}$, PGD-$L_{\infty}$ and Nominal training.
Results for the CIFAR data set are shown in Table 15. We observe that the
proposed methods ARUB-$L_{1}$ and ARUB-$L_{\infty}$ yield the highest
adversarial accuracies with respect to PGD-$L_{2}$ attacks. For the MNIST data
set and the FASHION MNIST data set (Table 16 and Table 17, respectively), we
see that ARUB-$L_{\infty}$ has the highest adversarial accuracies with respect
to PGD-$L_{2}$ attacks when $\rho\leq 0.1$, whereas PGD-$L_{\infty}$ does best
for larger values of $\rho$.

$\boldsymbol{\rho}=$ | 0.000 | 0.010 | 0.020 | 0.030 | 0.100 | 1.000 | 3.000 | 5.000 | 10.000
---|---|---|---|---|---|---|---|---|---
ARUB-$L_{1}$ | 71.17 | 70.98 | 70.70 | 70.51 | 69.88 | 60.31 | 44.69 | 35.35 | 17.89
ARUB-$L_{\infty}$ | 72.03 | 71.76 | 71.45 | 71.25 | 69.84 | 56.13 | 43.98 | 34.26 | 16.25
PGD-$L_{\infty}$ | 70.78 | 70.55 | 70.43 | 70.39 | 69.65 | 59.96 | 43.59 | 29.10 | 16.21
Nominal | 71.29 | 71.13 | 71.05 | 70.66 | 69.38 | 48.16 | 16.05 | 10.04 | 10.04
Table 15: Adversarial Accuracy for CIFAR with CNN architecture and PGD-$L_{2}$ attacks. $\boldsymbol{\rho}=$ | 0.000 | 0.010 | 0.020 | 0.030 | 0.100 | 1.000 | 3.000 | 5.000 | 10.000
---|---|---|---|---|---|---|---|---|---
ARUB-$L_{1}$ | 91.25 | 91.09 | 90.78 | 90.55 | 89.02 | 82.46 | 68.95 | 54.80 | 33.20
ARUB-$L_{\infty}$ | 91.37 | 91.13 | 91.05 | 90.86 | 90.59 | 82.81 | 71.45 | 60.23 | 30.86
PGD-$L_{\infty}$ | 90.59 | 90.55 | 90.43 | 90.23 | 89.30 | 84.45 | 73.20 | 67.54 | 57.77
Nominal | 91.02 | 90.90 | 90.62 | 90.43 | 89.02 | 77.85 | 56.72 | 40.51 | 10.16
Table 16: Adversarial Accuracy for Fashion MNIST with CNN architecture and PGD-$L_{2}$ attacks. $\boldsymbol{\rho}=$ | 0.000 | 0.010 | 0.020 | 0.030 | 0.100 | 1.000 | 3.000 | 5.000 | 10.000
---|---|---|---|---|---|---|---|---|---
ARUB-$L_{1}$ | 99.38 | 99.38 | 99.38 | 99.34 | 99.30 | 98.32 | 91.68 | 70.51 | 33.83
ARUB-$L_{\infty}$ | 99.38 | 99.38 | 99.38 | 99.38 | 99.30 | 98.79 | 95.70 | 87.46 | 47.93
PGD-$L_{\infty}$ | 99.22 | 99.22 | 99.14 | 99.14 | 99.02 | 98.55 | 96.37 | 91.29 | 55.39
Nominal | 99.30 | 99.26 | 99.26 | 99.26 | 99.18 | 97.93 | 89.69 | 60.00 | 9.34
Table 17: Adversarial Accuracy for MNIST with CNN architecture and PGD-$L_{2}$
attacks.