SmoothHess: ReLU Network Feature Interactions via Stein's Lemma
Authors: Max Torop
Aria Masoomi
Davin Hill
Kivanc Kose
Stratis Ioannidis
Jennifer Dy
Venue: NeurIPS 2023 poster (2023)
Keywords: Interpretability
Feature Interactions
Stein's Lemma
TL;DR: We use Stein's Lemma to compute a smoothed Hessian for interpreting feature interactions in ReLU networks.

Abstract: Several recent methods for interpretability model feature interactions by looking at the Hessian of a neural network. This poses a challenge for ReLU networks, which are piecewise-linear and thus have a zero Hessian almost everywhere. We propose SmoothHess, a method of estimating second-order interactions through Stein's Lemma. In particular, we estimate the Hessian of the network convolved with a Gaussian through an efficient sampling algorithm, requiring only network gradient calls. SmoothHess is applied post-hoc, requires no modifications to the ReLU network architecture, and the extent of smoothing can be controlled explicitly. We provide a non-asymptotic bound on the sample complexity of our estimation procedure. We validate the superior ability of SmoothHess to capture interactions on benchmark datasets and a real-world medical spirometry dataset.

SmoothHess: ReLU Network Feature Interactions via Stein's Lemma
Max Torop," ' Aria Masoomi," 1+ Davin Hill,' Kivane Kose, 2 Stratis Ioannidis, Jennifer Dy' Northeastern University, 2 Memorial Sloan Kettering Cancer Center (torop.m, m, masoomi, @/@oothheateern edu, {dhill, ioannidis, jdy) neu. edu, (kosek)mmskcc.oo org
Abstract
Several recent methods for interpretability model feature interactions by looking at the Hessian of neural network. This poses challenge for ReLU networks, which are piecewise-linear and thus have a zero Hessian almost everywhere. We propose SmoothHess, method of estimating second-order interactions through Stein's s Lemma. In particular, we estimate the Hessian of the network convolved with a Gaussian through efficient sampling algorithm, requiring only network gradient calls. SmootHHess is applied post-hoc, requires no modifications to the ReLU network architecture, and the extent of smoothing can be controlled explicitly. We provide non -symptotii bound on the sample complexity of our estimation procedure. We validate the superior ability of SmoothHess capture interactions on benchmark datasets and real-world medical spirometry dataset.
1 Introduction
As machine learning models are increasingly relied upon in a variety of high-stakes applications such credit lending [50, 42] medicine [10, 67, 15], or law [38], it is important that users able to interpret model predictions. To this end, many methods have been developed to assess the importance of individual input features in effecting model output [53, 76, 81, 66, 73]. However, one may achieve a deeper understanding of model behavior by quantifying how features interact to affect predictions. While diverse notions of feature interaction have been proposed [53, 57, 56, 63, 25, 77, 84], this work, we focus on the intuitive
definition of feature interaction [25, 4, 47, 85, 18, 39].
Given function and point, the Interaction Effect [47] between given set of features on the output the partial derivative the function output with respect to the features; intuitively, it represents the infinitesimal change in the function engendered by joint infinitesimal change in each chosen feature. The Interaction Effect derives in part from Friedman and Popescu [25], who define the global interaction between set of features over the data distribution as the expected square of the partial-derivatvee with respect to those features. As in prior works [47,39] we eschew the expectation focus on local interactions occurring around given point x and avoid squaring partial-derivatives to maintain the directionality of the interaction. We focus on pair-wise feature interactions, which, in the context of the Interaction Effect,, are the elements of the Hessian.
ReLU networks are a popular class neural networks that use ReLU activation functions [28]. The of ReLU has desirable properties, such the mitigation of the vanishing gradient problem [28]. and the sole activation function used in popular neural network families such as ResNet [32] and VGG [72]. However, ReLU networks are piece-wist linear [58] and thus have everywhere (a.e.), posing problem for quantifying interactions (see Figure (aa)..
Hessian almost
"Equal contribution
37th Conference on Neural Information Processing Systems (NcurIPS 2023).

---

x2)
(a)
(b)
5
Figure 1: (a) An exemplar illustration of simple 5-hidden-layee ReLU network f R - R. Note that / is piece-wiss linear and thus has (1x(x) 0 a.c. (b) The ReLU network convolved with qu.3z B + R, the density function (V 0.37). This function longer piece-wise linear and admits non-zero higher order derivatives.
A common approach for estimating the Hessian in ReLU networks is to take the Hessian of smooth surrogate network which approximates the original: each ReLU is replaced with SoftPlus, smooth approximation to ReLU [28, 22]. before differentiating [39, 20J. However, this approach affords one only coarse control over the smoothing as each internal neuron is smoothed, leading to unwieldy asymmetric effects, can be seen in Figure 2.
In this work, we propose SmoothHess: the Hessian of the ReLU network convolved with a Gaussian. Such a function is more flexible smooth surrogate than SoftPlus as smoothing is done on the network outpur, and the covariance of the Gaussian allows one mediate the contributions of points based on their Mahalanobis distance. Unfortunately, obtaining the Hessian of the convolved ReLU network is impossible using naive Monte-Carlo (MC) averaging. However, by proving an extension Stein's Lemma [78, 49], we show that such a quantity can be efficiently estimated using only gradient calls on the original network. For an illustration of ReLU network smoothing, see Figure 1(b).
Our main contributions are as follows:
We propose SmoothHess, the Hessian of the Gaussian smoothing a neural network, as a model of the second-org -oceed-order feature interactions.
We derive a variant of Stein's Lemma, which allows one to estimate SmoothHess for ReLU networks using only gradient calls.
We prove non-asymptotic sample complexity bounds for our SmoothHess estimator.
We empirically validate the superior flexibility of SmoothHess to capture interactions on MNIST,
FMNIST, and CIFARIO. 10. We utilize SmoothHess to derive insights into a network trained on real-world medical spirometry dataset. Our code is publicly available."
The remainder of this paper is organized as follows: In Sec. 2, we summarize the gradient-based methods for feature importance and interactions that are most related to our work. In Sec. 3, we provide a technical preliminary covering the definitions and techniques underlying both our method and competing works. Next, in Sec. 4, we introduce estimate it, and provide sample complexity bounds for our estimation procedure. In Sec. 5, we experimentally demonstrate the ability of SmoothHess to model interactions. Finally, in Sec. summarize our method and results, followed by discussion of limitations and possible solutions.
method, SmoothHess, explain how to
Related Work
Feature Importance and First-Order Methods: Methods that quantify feature importance fall into
two categories: (i) perturbation-based methods (e.g.. [53, 66, 17]), which evaluate the change in model outputs with respect to perturbed inputs, and (ii) gradient-based methods (e.g., [72, 76, 81].. which leverage the natural interpretation of the gradient as infinitesimally local importance for given sample. Most relevant our work are gradient-based approaches. The saliency map, as defined in [72], is simply the gradient of model output with respect to the input. Several variants are developed address the shortcomings of the saliency maps. SmoothGrad [76] was developed
"https://githuu htt.ha/ww..cooaaaoooot/tttttttee

---

SH: vf(f*
(,) (F(f*,n)0);
log o²
SP Hess V2f((0)
log 8
(a) 10
(b)
LOING a 00 00
Smoother
1000
Smoother
Vit a2100,
2 1o
10g
Figure 2: Estimated Hessian element between features and 2 f: + R trained to memorize the Four Quadrant toy dataset. (a) SmoothHess (SH) is estimated with isotropic covariance If1 using granularly sampled o² E (Ie-3..... 1}. Aside from minute o² 1e-2.5, where hyper-local noisy behavior is captured, we have (f (zo) (z=))) = (5 15++++ + the memorized 'ground truth" off off-diagonal Hessian element over the four quadrants. This indicates symmetry in the weighting of the contributions from points around at every level of smoothing. (b) The Hessian of the SoftPlus smoothed function fa (SP Hess) computed using granularly sampled average value 2.5 is not achieved any value 3, aside from briefly between log1o 3 and log $ indicating that SoftPlus fails to incorporate the information around xo in a symmetric manner at every level of smoothing.
(0, 0)" for 6-layer ReLU Network
10)/4 = 2.5, the average of
fle-1,.... le4]. The
to address noise by averaging saliency maps (see also See. 3), and comes with sample complexity guarantees [3]. Sundararajan et al. [81] introduce Integrated Gradients, the path integral between an input and an uninformative baseline, This is extended to the Shapley framework by Erion et al. [24]. The Grad-CAM line of work [93, 70, 61] is similar in nature to the methods above, with the key distinction that importance is modeled over internal (hidden) layers.
Feature Interactions: A variety of methods have been proposed to estimate higher-order feature interactions, which can again be separated into perturbation-based [54, 82, 86, 57, 92] and gradient- based approaches. Among the latter, Tsang et al. [85] propose Gradient-NID, which estimates feature interaction strength the corresponding Hessian element squared. Janizek al. [39] propose Integrated Hessian, which extends Integrated Gradients to use path-integrated Hessian. Lerman et al. [47] propose Taylor-CAM, higher-order generalization global interactions for Bayesian neural networks in terms of the expected value of the Hessian over the data distribution. For classification ReLU networks, the CASO and CAFO explanation vectors exploit feature interactions in the loss function using Hessian estimates [74].
Grad-Cam [70]. Cui et al. [18] quantify
Unfortunately, due their piecewise linearity, existing fessian-based interaction methods cannot be readily applied to ReLU networks. Janizek et al. [39] replace each ReLU activation with SoftPlus post-hoc, before applying Integrated Hessians. Similarly, Tsang et al. [85] apply their method to networks with the SoftPlus activation instead of ReL.U. For regression tasks, Lerman et al. [47] replace ReL.U with the smooth activation function GELU [34] before training. Although SoftMax outputs and the cross-entropy loss admit higher order derivatives, pre- post-hoc smoothing, as above, is necessary for finding interactions affecting logits, internal neurons, or regression outputs. Indeed, while Singla et al. [74] estimate interactions the original ReLU network, they are only with respect to the loss function.
In contrast, we propose method for quantifying feature interactions that works with any ReLU network post-hoe without requiring retraining modifications the network architecture. also can be directly estimated with respect to model as well as intermediate layer outputs. Furthermore, our experiments in Sec. 5 show the superior ability of our method to model logit outputs, internal neurons, and SoftMax probabilities as compared to SoftPlus smoothed network.
3 Technical Preliminary
We denote ReLU network by F Rd R where c= in the case of regression. We denote the function which we wish to explain as f Rd - R, an arbitrary neuron f!0) (the hth neuron in the 1th layer) in our ReLU network, a SoftMax Probability
ReLU Network Background:
3

---

for some class & E cl). We denote sample space X' C R and point Tu E X' for which we wish to capture feature interactions affecting f (xo). A general L-hidden-layer ReLU Network F = f((+)) Rd + R may be defined recursively by [33, 58]:
(Fo) = /( (za) + 80) for 1 = 1 1....... I. 1,
(1a)
( max(0, (xo)) for ( = I...L,
(lb)
with (x0) = 2g E R. For a given layer the dimension (number of neurons) is defined as ny E N, Le. f (xo), g0) (20) E R (with nop = d) with weight and bias W E & b()) E R" As
stated above, ReLU networks are piecewise linear [58]. Specifically, each neuron f. R + R, E
(1....,++1],. E {1..... m. is a piecewise linear function, corresponding to a finite set of K E N convex polytopes = (Q.) Q S IRd which form a partition of Rd [30, 31, 58, 33, 71, 64, 7]. See Figure ((a) for an example of ReLU network.
SmoothGrad [76] is an extension of the saliency map (x.( [72] developed to reduce noise. SmoothGrad is average over saliency maps; formally, it is Ess V. (zo + d)] where 8 ~ (NO,E) Σ) and the covariance Σ E Raxd is a hyperparameter usually chosen to be isotropic 0² I²I,o Ø > O. SmoothGrad is estimated by sampling n perturbation vectors {6, 16. 6, NO., Vi E (1.....n) n) and averaging via GSG (Zo, f, Σ) = 1 n Ex= ( + 6.). Unfortunately, (z00 s Es I 6)] is not useful for quantifying
SmoothGrad:
the analogous Hessian average Ei= V feature interactions in ReLU networks as I& (V + 5)] = O.
SoftPlus: A common approach to assessing higher-order derivatives for ReLU networks is by
differentiating a smooth surrogate fo R - I where a SoftPlus replaces every ReLU in f [28, 39, 20, 21J. The SoftPlus function s (0 = 3 log + e B > 0, is a smooth approximation of ReLU, where 3 is a parameter inversely proportional to the level of smoothing [28, 22]. SoftPlus approaches ReLU and Sa and V/s approximate ʃ and V with arbitrary accuracy as 8 - oo [20]. Dombrowski et al. [20] establish empirical and theoretical connections between V= and SmoothGrad. They observe visual similarities between va (2ro) and SmoothGrad at appropriate values. Further, in the simple case of single neuron with no bias, they prove that replacing ReLU with SoftPlus corresponds to a Gaussian- smoothing of the gradient. However, this theory does not extend to networks with multiple neurons and layers.
Stein's 's Lemma: Stein's Lemma [78] is central to our analysis. We relate variant of Stein's Lemma [51], which extends the results of Stein [78] to hold for multivariate normal distributions with arbitrary covariance matrices:
Lemma 1. (Stein's 's Lemma (511) Given xo E Rd, covariance matrix E E Raxd multivariate normal random vector d E Id distributed from ~ (0, Σ) and almost everywhere differentiable functioon g Rd + for which E3Il[Va9(a) + {)] < oc for each E {l..... d}, then
Eg - 5g( + 8)] = Ev[V g0z + 5)].
(2)
A zero-th order oracle associated with a function Rd - IR is an oracle which, when provided with any input i E Rd, returns g(x). Likewise, first order oracle (or gradient oracle) returns Vzg(z). In the context of this work, where is ReLU network, zero-th order and first order oracle calls amount to network forward passes and backpropagation calls, respectively.
Given one has access to as a zero-th order oracle, the LHS of Stein's 's Lemma may be estimated by sampling set of perturbations (o.) 8, ~ (O, E), querying the zero-th order oracle g(zo + 6:) for each o, and MC-estimating, i.e:
Go 9. 2)
x x
16g(20 + 5.) E3 ōg(20 + 8)] Lemma
E V20 + +)].
(3)
n
Such an approach is useful for estimating E. VV.g(ar, (2n + 6)] the RHS of Eq. (3) when the gradient of y is impossible or expensive to obtain but g itself may be efficiently queried as zero-th order oracle (see e.g., [44, 48, 60)). Extending Liu [5]]'s work, Lin al. [49] present second-order
variant of Stein's Lemma expressing the expected Hessian in terms of the gradient:

---

Lemma 2. (First-Orer Ordsrr Oracle Stein's Lemma [49]) Given 2o E R covariance matrix E E Bdxd multivariate normal random vector å E R distributed from d ~ (0, Σ) and continuously differentiable function g(2) IR IR with locally Lipschitz" gradients Vg R² + R2 then
E8 Σ V[V 200000000 + 5)]² = E& E 6)]-
Complexity bounds have been derived for similar identities, which express the Hessian using zero-th order information [6, 94, 23]. However, Lemma 2 faily for ReLU networks. This is precisely because ReLU networks are piecewise linear and, therefore, are not continuously differentiable. In the next section, we directly address this through our method for estimating smoothed Hessian for ReLU networks.
4 SmoothHess
Our main contributions are: (1) We propose SmoothHess, the Hessian of the network convolved with Gaussian, for modeling feature interactions. (2) We use Stein's Lemma prove that SmoothHess may be estimated for ReLU networks using only gradient oracle calls. (3) We prove non-asymptotic sample complexity bounds for SmoothHess estimator.
Gaussian Convolution as a Smooth Surrogate: An alternative smooth-surrogate to Sa is h
Rd TR, the convolution of f with Gaussian:
h = (S * (20)) = - f()) zo)dz,
"EER"
where Rddd is a covariance matrix, qz(2 20) (2 exp( exp(-3d(<oo, (0(o 2)x) is the density function of the Gaussian distribution E), is the determinant and d(ro, 2) (xo- --) (x0- - 2) E the Mahalanobis distance between Fg and
The Gaussian Gaussian-smoothed function h is infinitely differentiable and does not suffer from the limita- tions of surrogates obtained from internal smoothing. Here, smoothing is done the output image of and thus the relationship between h (xto) f(=). = and xo is made explicit by Eq. (5): the relative the exponentiated negative half Mahalanobis distance } d z)x. The ability adjust Σ gives user fine-grad fine-grained and localized control over
contribution of f=) (20) is proportional
smoothing: the eigenvectors and eigenvalues of Σ respectively encode directions of input space and corresponding locality for their contribution to hy
We define SmootHHess as the Hessian of hy
Definition 1. (SmoothHess) Given ReLU network J :Rd TR, point to explain E IRd,, covariance matrix E Raddd and qr. Rdd R. the density function of Gaussian distribution N N(0, E), then SmoothHess is defined to be the Hessian of convolved with qr evaluated at zo:
-f SSER f(2) (000)((2=))dd
Vahy,(0z)) = vfllf { (=)
(6)
Well-known properties of the Gaussian distribution may be used to encode desiderata into the convolved function and, accordingly,
SmoothHess through the choice of the covariance. For instance, it is known that as d - oo the isotropic Gaussian distribution N/ (f, converges to U(Sd)1). a uniform distribution over the radius a d sphere [88]. Thus, given large enough d, as is commonly encountered in deep learning datasets, one may choose Σ (r/vd)/ to approximately ensure that ( incorporates information from the radius r sphere around Øg. We exploit and validate this intuition in our experiments (see Table 1).
Finally, we must highlight the strong connection between SmoothHess and SmoothGrad [76], which Wang et al. [89] prove is equivalent to the gradient of the same smooth surrogate: V v, Thus,
most general form, Lemma 2 holds for functions which are continuously differentiable and have locally ACL gradients. Locally ACL functions are functions which are absolutely continuous on almost every straight line, mild condition which is satisfied by both locally Lipschitz and continuously differentiable functions [49, 46, 68). As ReLU networks are locally Lipschitz [29], they are locally ACL.
5

---

SmoothGrad and SmoothHess together define second-order Taylor expansion of f, at zo. which can be used second-order model of f around Zg.
SmoothHess Computation via Stein's Lemma: We relate our method for estimating SmoothHess for ReLU networks. As stated above, Lemma 2 does not hold for ReLU networks. However, extend the arguments of Wang et al. [89] and Lin et al. [49] to show that the LHS from Lemma 2 is equivalent to SmoothHess for all Lipschitz continuous functions":
Proposition 1. Given E TRd,, I.-Lipschitz continuous function 9 R - RR covariance matrix E E
Rdxd and random vector E IR distributed from ó ~ (O, 2)) with density function qr R + IR then
Es - 8 Vag(o + +8)7 = 22 qx)) = (ihgo(00))
(7)
where: denotes convolution.
In other words, even though Lemma 2 does not hold for ReLU networks, Stein's Lemma is indeed evaluating Hessian: namely, SmoothHess given by Eq. (6). The proof consists of moving the Hessian operator on the RHS of Eq. (7) into the integral defined by (g * (zn). The resulting expression is simplified into form for which Lin et al. [49] prove is equivalent the LHS of Eq. (7). The proof is provided in App. B.
Proposition 1 opens up the possibility of MC-estimate of SmoothHess that only require access to a first-order oracle f. This is computed by sampling a set of n E N perturbations {6 5, ~ (0, E), and for each 6, querying V + 6,) before taking the outer product 6, (x + 6,))² Rdxd Monte-Carlo averaging and finally symmetrizing:
(oo, f.E) = E
5, V z + o,)"",
(8a)
An f,2)= (1 2. + fr." (xo, f.).
(8b)
Each first-order oracle call has the same time complexity one forward pass and may be computed efficiently in batches using standard deep learning frameworks [62, 1, 13]. As is finite, (20, f. E) is not guaranteed to be symmetric in practice. Thus, we symmetrize our estimator in Eq. (8b). A straightforward consequence of Proposition 1 is that lim (2o. f. Σ) = (Zo)) which we formally show in the Proof of Theorem 1 in App. C.
Estimation of SmoothGrad may be amortized with SmoothHess, obtaining V y (20) at a signifi- cantly reduced cost. The main computational expense when estimating SmoothGrad is the querying of the first-order oracle. As this querying is part of SmoothHess estimation, the gradients which we compute may be averaged at minimal additional cost of O(nd) to obtain SmoothGrad estimate. Likewise, SmoothHess may be obtained at a reduced cost of O (ndd),, ), the cost of the outer products, during SrodthGrad estimation. For details of our algorithm, App. D.
Last, we prove non- asymptttiiiic bounds for our SmoothHess estimator:
Theorem 1. Let f R R. be a piece-wise linear function over a finite partition of Rdd Let Fp E Rdd and denote { a set of n iidd random vectors in Rdd distributed
from 6, ~ NNO,E) Given Hn Σ) as in Eq. (8), for any fixed E,y E (0, 1]. given n 2 $ [max V oF log $2, (C- Vd +
log 1)2)] then
P -1|
H112 >0) S7
where H = V (zo)), C+,00 C,, 0 are constants depending on the function f and covariance 2 and qr : Rd + R is the density function N0, 2).
We elect to the non-asymptotic bounds presented in Theorem 3.39 of Vershynin [87], which hold for sums of outer products of a sub-gaussian random vector with itself. This poses challenge as H, is the sum of outer products of two sub- oub-gaussian random vectors which are not necessarily equal. To deal with this issue, we separately prove asymptoti bounds for the outer products
*All ReLU network outputs, internal neurons, and SoftMax probabilities are Lipschitz continuous [29, 26].

---

Dataser
MNIST
FMNIST
CIFARIO
Function
Class I Logit (1)
Int. Neuron
Class
Logit
Int. Neuron
Class Logir (1)
Int. Neuron (8)
0.25 0.50 1.00
0.25 0.50 1.00 0.25 0.50 1.00
0.25 0.50 1.00 0.25 0.50 1.00
0.25 0.50 1.00
SH++G (Us)
9.6-77 7.8-6 6.7e-5
4.9e-8 4.0e-7 3.3e-6 6.5e-7 4.0e-6 4.3e-5
2.0c-8 18c-7
9.8e-4 2.2e-2
8.1e.7
16e-4
[76]
4.5e.6 4.1e:5
2.1e.7
1.7e.6 1.5e-5 3.0e
2.7e:5
1.0e77 9.0e:7 7.0e-6
1.3e:2
8.6e
1.3e:5 1.le.4 8.3e
+ G)
1.2c-6 9.66-6
8.16-5
5.5e-8
4.4c-7 3.70-6 9.6e-7
7.5c-6 6.5o-5
3.0c-8
2.1e-7 1.8c6)
2.1c-3
3.3-2 2.5c-1
1.1e-5 1.0c-4 7.004
SP
4.6e-6 1.1e-64.1e-5 3.9e-
2.1e-7
1.7c-6 1.50-5 3.2e-6 2.8c-5 2.6c-4
1.0c-7
8.5e-7 7.20-6
1.3e-2 9.0e-2 5.2e-1
5.1e-5
2.9c-4 1.60-3
G
4.2e-3 1.7e-2 6.7e-2
2.0e-3 7.0e-3 2.9e-2 3.8e-3 1.5e-2 6.0e-2 1.0e 4.0ee4 1.8e-3 3.0e-1 1.2e-0 5.0e
5.00.099.0
3.5e-3 3.55-
1.4e2
Table I: Average PArsE at three radii Results are calculated for SmoothHess + SmoothGrad ISH+SG,Us) SmoothGrad (SG) SoftPlus Gradient (SP G) SoftPlus Hessian + SoftPlus Gradient (SP (H + G)) and the vanilla Gradient (G). Results are provided for the predicted class logit, and the penultimate neuron maximally activated by the "three." dress and cat classes for MNIST. FMNIST and CIFARIÓ respectively. While SP (H + G) and SP G results are reported using the best / chosen from over 100 values based on validation set performance, only 3 values were checked for SH SG and SG based upon a = <) Vdyd. SH SG outperforms the competing methods for all 18 permutations of dataset, function and e.
of the positive and negative eigenvectors of the summands. This is accomplished using an identity expressing the eigenvalues of the summands in terms of E-'8 and V.f(rr + 8), which we state and prove in App. C. The proof is completed by applying the triangle inequality and union bound to combine the two bounds into Eq. (9). Our proof is included in App. C.
5 Experiments
5.1 Experimental Setup
Experiments were conducted on real-world spirometry regression dataset, three image datasets (MNIST [45], FMNIST [90] and CIFARIO [43]), and one synthetic dataset (Four Quadrant). The Four Quadrant dataset consists of points x' E R²² sampled uniformly from the grid 2, x [2, 2] C TR22, with spacing of 0.008. The label of point y(x') = K E IR is chosen based upon its quadrant, with K 5, 3, 12, - 10 for Quadrants 1, 2, 3 and 4 respectively. We train a 5-layer network on MNIST and FMNIST a ResNet-18 [32] on CIFARIO 10 and a 6-layer network on Four Quadrant. A 10-layer 1-D convolutional neural network was trained
Datasets and Models:
the spirometry regression dataset. Additional model, hyperparameter and dataset details are outlined in App. E.
Hardware: Experiments were performed on an internal cluster using NVIDIA A100 GPUs and AMD EPYC223 7302 16-Core processors.
Metrics: Gradient-based explainers are evaluated using a Taylor expansion-like function as a proxy. We estimate a gradient Hessian pair G( E Rdd H (( E Rddd for some function : Rdd -> IR around point To. Here g is either smooth- moothh-uurooatee of f or f itself. We define the following Taylor expansion-like function:
1 (A 2o) 20) -10)-
(A) (,, G mo) + 2
(10)
For readability, we remove the dependence on To from ((2..)) and H (mo.a)- Eq. (10) is almost equivalent to the Taylor-expansion of g Fo, with the key difference that the zero-th order term is ((0) as opposed g(zr)) This is done isolate the impact of the explainers Gg and Hl Setting H = yields J.o.(G,.0(4) (A)) first-order Taylor expansion-like function. For brevity, we refer to F f H,, simply as Taylor expansion below. We use the following two metrics to quantify the efficacy of SmoothHess:
We introduce the Perturbation dean-Squared Error (PMSE) to assess the ability of il and/or G to capture the behaviour of f over given neighborhood. Given point xo and radius 8 > 0, the P PASE directly measures the fidelity of f G,, to when restricted to the ball B, (2o):
Perturbation MSE:
[ (2EB.
PMSE(50, G,, (ol(Ba(zo)
(x') (, {(x)))
(11)
A low PMSE value indicates that R, is good fit f when restricted to (zo),, and thus that Go and Hg capture the behaviour of of over (zo))). PMSE may be estimated for point
7

---

Dataset
MNIST
FMNIST
CIFARIO
Attack Mognizude
0.25
0.50
0.75
1.25
1.75
0.25
0.50
0.75
1.25
1.75
0.1
0.2
0.3
1.0
SH+SG (Us)
93.0
80.3
48.0
10.5
2.0
79.5
46.8
25.0
3.5
0.0
62.5
38.5
26.5
15.0
4.5
[76]
81.8
48.8
79.5
49.3
26.3
27.5
0.0
SP (II+G)
81.8
51.5
79.8
51.0
27.5
64.5
31.0
7.5
SP
93.3
82.3
53.8
16.3
5.0
79.8
51.5
29.5
7.8
1,0
66.5
47.5
36.0
29.5
8.5
[73]
93.3
82.8
56.0
18.5
8.8
80.3
52.3
31.8
11.0
2.5
69.0
51.5
41.0
34.0
21.5
Random
99.8
99.5
99.0
99.0
98.8
99.3
98.0
97.3
95.5
93.8
100.0
99.5
99.0
98.5
96.5
Table 2: Post-hoc accuracy of adversarial attacks performed on the predicted SoftMax probability, at five attack magnitudes Lower is better. Results for SmoothHess SmoothGrad (SH SG, Ours), SmoothGrad (SG), SoftPlus Gradient (SP G) and SoftPlus Hessian + SoftPlus Gradient (SP (H + G)) are reported using parameters 0²1 E Raxd and / > O chosen based upon performance on a held-out validation set. We additionally compare with the vamilla (unsmoothed) Gradient (G). First order attack vectors are constructed by scaling the normalized gradient by and subtracting from the input. Second order attack vectors are found by minimizing the corresponding second-order Taylor expansions.
To by sampling a set of n E N points {= uniformly from Bg computing the errors and MC-averaging: 1 D (x() - f(x1))² = PasSE(Fo. f. f.f.fooooe G8,
Adversarial Attacks: Given & and H, one can use Fz., G. to construct adversarial attacks of any desired magnitude > O. We denote A E Rd as the perturbed input, and corresponding attack vectors by A zo. Given only Gg (i.e. Hg = 0), the first-order Taylor expansion attack yields: A Io G5 =C,,\Il., /|1G Otherwise, the attack may be framed as the solution quadratic optimization minimizing the second-order Taylor expansion:
min AERd G20.G,, (4), s.t. | - - S..
(12)
Although this optimization is non-convex (as Hg is not guaranteed to be positive-semi-defnie)) can he solved exactly [12]. For implementation details, and a discussion of the similarities and differences with Singla al. [74], see App. D. We set g to be the smoothed version of the predicted class SoftMax
probability function. Given set of test points { 11 E N. we validate the efficacy of our attacks using the average post-hoc accuracy metric n E I[argmax, Fxx;) = argmax, F(A;)), where I denotes the indicator function, and A the output after z; is attacked.
Setup: We provide the details for our experiments below:
Four Quadrant: We train the network to memorize the Four Quadrant dataset. We measure the interactions at æ0 (0, 0)* using both SmoothHess and SoftPlus Hessian, estimated with granularly sampled d² E {le-3..... 1} and /3 E {le-1.... 4}.
P PMSE:
We set g to be the smoothed version of a predicted class logit or penultimate neuron. The penultimate neuron is chosen to be the maximally activated neuron on average over the train data by the "three." dress and cat classes for MNIST, FMNIST and CIFARIO, respectively. We compute PMSE for three neighborhood sizes E (0.25, 0.50, 1.0). While over one-hundred values of 3 are checked on validation set before selection, only three values of and SmoothGrad, based the common-sense criterion = e/vd: E {e/2 3e/4 e e/v outlined in Sec 4. The standard deviation of PASE results is reported in App.
checked for SmoothHess
Adversarial Attacks: Adversarial attacks are performed after selecting the best B and from held-out validation set. We refrain from using the criterion = e/ d for adversarial attacks, as they rely upon the extremal, as opposed to average, behavior of f. We check between ss 10 and 30 values of a and 3 the held-out validation set before selecting the values resulting in the most effective attacks. Attacks are performed using magnitudes= e E {0.25, 0.50, O. 75, 1.25, 1.75) for MNIST and FMNIST and {0.1, 0.2, 0.3, 0.4, 1.0} for CIFARIO, which is easier to successfully attack at lower magnitudes due to it's complexity. For more details see App. D.
Competing Methods: We compare SmoothHess with other gradient- brased methods that model f locally; i.e. those which can associated with the Taylor expansion around a smooth surrogate or f itself. Specifically, we compare with the Gradient and Hessian of SoftPlus smoothed network Ss [20, 39], SmoothGrad [76] and the vanilla (unsmoothed) Gradient. We also compare adversarial attacks with random vectors scaled to the attack magnitude as baseline.
8

---

Sample 1: Early Cough
Sample 2: Late Cough
SmpathGrad
SmpothHless
SmoathGrad
SmoothHess
0.04
0.010
0.04
3,010
4000
3000
Interactinn
4002
0.003
Votume 1P2
4000
5.005
Smood Hose
2000
2000
2000
3,000
1000
0.520
1000
0.005
2000
-2.02
1000
0.005
6
0.04
-.010
0.030
Time
Tims Sal ometryy Curve
5G
Time [4]
Time (1)
Figure 3: Evaluation of network predictions for two spirometry samples (Left and Right), using SmoothGrad and SmoothHess. Spirometry curves are plotted in blue. SmoothGrad (green) and SmoothHess (red) attributions are grouped into 0,5s intervals and plotted as bars for cach interval. The plotted SmoothHess values are interactions with respect to the first 0.5s time interval. The red arrows point to the features that exhibit the largest interactions (endpoint) with the initial 0.5s (source). Sample exhibits coughing within the first 4 seconds, as indicated by plateauing spirometry curve. We observe that the strongest interactions for Sample occur this initial cough-induced plateau. In contrast, Sample 2 exhibits no coughing within the first 4 seconds; the strongest interactions occur later near small fluctuations in the spirometry curve.
Additional experiments are provided in App. F. We compare the efficacy of SmoothHess with the unsmoothed Hessian (for adversarial attacks on SoftMax) and Swish [65] smoothed networks, both of which our method outperforms. We also run the PASE experiment using a ResNet101 model, validating the superior ability of SmoothHess to capture interactions in larger networks.
5.2 Results
Symmetric Smoothing- Four Quadrant Dataset: We investigate the ability of both SmoothHess and the SoftPlus Hessian to capture local feature interactions symmetrically. Results are shown in Figure 2. At each value of o2, aside from extremely small 6² < 1e-2.5, the off-diagonal element of SmoothHess is approximately equal to 2.5, the average interaction over the quadrants. The off- diagonal element of the SoftPlus Hessian is essentially never equal to 2.5. The results for SmoothHess follow from the fact that an isotropic covariance was used, the rotation invariance of which weights points that are equidistant from x equally. The SmoothHess off-diagonal element is not 2.5 at small a because, despite the network being trained to memorize the data, such small neighborhood will inherently reflect noise. The inability of SoftPlus to capture interactions symmetrically follows from the fact that smoothing is done internally on each neuron, which does not guarantee symmetry.
Perturbation Mean-Squared-Erron Error: We show SmoothHess can be used capture interactions occurring in a variety of neighborhoods around x. It is shown in Table 1 that SmoothHes, + SmoothGrad (SH SG) achieves the lowest PMSE for all 18 combinations of dataset, function and neighborhood size e. We emphasize that this is despite the fact that only three values of a were validated for SG + SH and SG on a held-out set compared over meehundredd values of for SP (H G) and SP G. This indicates the superior ability of SmoothHess to model f over a diverse set of neighborhood sizes e. Further, it can be seen that second-order methods SH + SG and SP (H + G) achieve significantly lower P than their first-order counterparts, SG and SP G, respectively. sometimes by an order of magnitude. This confirms the intuition that higher order Taylor expansions of the smooth surrogate provide more accurate models of f. Interestingly. while SH + SG is clearly superior to SP (H + G), we see that SG and SP G are tied in many cases. This could indicate that the symmetric properties of Gaussian smoothing are comparatively more important when higher-order derivatives are considered. However, more investigation is needed.
We use the interactions found by SmoothHess to generate adversarial attacks on the classifier. We see from Table 2 that the Taylor expansion associated with SmoothHess and SmoothGrad generates the most powerful adversarial attacks for all datasets and attack magnitudes E. aside from CIFARIO at 1.0 where both methods successfully change most class predictions. The superiority of the SmoothHess attacks indicates that the Gaussian smoothing is capturing the extremal behavior of the predicted SoftMax probability more effectively than the other methods. One hypothesis for SmoothGrad outperforming SmoothHess in the largest neighborhood for CIFARI0 that the network behavior is highly complex over a large area, and adding higher order terms decreases performance.
Ad Adversarial Attacks:

---

Qualitative Analysis of FEV, Prediction using Rejected Spirometry:
spirometry test is common procedure used evaluate pulmonary function, and is the main diagnostic tool for lung diseases such as Chronic Obstructive Pulmonary Disease [40]. During an exam, the patient blows into a spirometer, which records the exhalation as volume-time curve. Recent works have investigated the use of deep learning on spirograms for tasks such as subtyping [11], genetics [16, 91]. and mortality prediction [36]. Additional background on spirometry is outlined in App. E.
Traditionally, exhalations interrupted by coughing are discarded for quality control reasons. We train a
Kernel Width: 200
Kernel Width: 20
ant b0st.
3075
CNN raw spirograms from the UK Biobank [79] to predict the patient's "Forced Expiratory Volume in I Second" (FEV 1), metric frequently used evaluate lung health, using efforts that were rejected due coughing.
<330
Fig. 3, we apply SmoothHess and SmoothGrad on two spirometry samples to under- stand their FEV1 predictions. In Sample I, coughing occurs within the first 4 seconds of the curve, as ev- idenced by the early plateauing of the curve which indicates pause in exhalation. In contrast, Sample exhibits no detectable coughing within the first 4 seconds. To improve the interpretability of the results, we group the features into 0.5 second time intervals. FEV, is traditionally measured in the initial seconds of the non-rejected samples (see App. E), therefore we calculate SmoothHess interactions with respect to the first 0.5 second time interval.
360
022RD1E SE SRFR
Figure 4: Heatmap of the SeoothHess interactions for the spirometry sample Fig. 3, calculated on CNNs with varying convolution kernel width.
FEV; is known be strongly affected by the presence of coughing [55]. The time intervals where coughing occurs, indicated by plateaus in the spirometry curves in Figure 3, should be an important signal for any model trained to predict FEV1. Indeed, we observe that for Sample I, the SmoothGrad attribution for the first 0.5s interval is relatively low, with strong interactions occurring at the cough- induced plateaus. In contrast, the first 0.5s interval for Sample 2 shows high importance, with lower magnitude interactions that may be indicative of small fluctuations in the volume-time curve.
In Figure 4 we present the SmootHHess matrix for Sample 1 in Fig. 3, applied on two CNN models with different convolution kernel width. Interestingly, the smaller kernel width constrains the interactions to features that are spatially close together. In contrast, the features in the large kernel model have long-nanging interactions. These results agree with the intuition that smaller kernel widths, which constrain the receptive field of neurons up until the final layers, may correspondingly limit feature interactions.
6 Conclusion and Future Work
We introduce SmoothHess, the Hessian of the network convolved with a Gaussian, as a method for quantifying second-order feature interactions. SmootHHess estimation, which relies only on gradient oracle calls, and cannot be performed with naive MC-averaging, is made possible by an extension of Stein's 's Lemma that we derive. We provide non-asymptotic complexity bounds for our estimation procedure. In contrast to previous works, our method can be run post-hoc, does not require architecture changes, affords the user localized flexibility over the weighting of input space, and can be run on any network output. We experimentally validate the superior ability of SmoothHess to capture interactions over a variety of neighborhoods compared to competing methods. Last, SmoothHess to glean insight into a network trained on real-world spirometry data.
The outer product computations in SmoothHess estimation have space and time complexity (dd)). When d > 0 this becomes expensive: for instance ImageNet [19] typically has d2 ~ 10000 This is common problem for all interaction effect methods. We leave as future work to explore computationally efficient alternatives or approximations to these outer products; a potential remedy is devise an appropriate power method to estimate the top eigenvectors of SmoothHess, rather than its entirety, which fits well our estimation via sampling low-rank matrices.
Limitations:
10

---

Acknowledgements
This project was supported by NIH grants ROICA240771 and U24CA264369 from NCI, in part by MSKCC's Cancer Center core support NIH grant P30CA008748 from NCI, and NIH 2T32HL.007227 41.
References
[1] Martín Abadi, Paul Barham, Jianmin Chen, Zhifeng Chen, Andy Davis, Jeffrey Dean, Matthieu Devin, Sanjay Ghemawat, Geoffrey Irving, Michael Isard, al. Tensorflow: system for large-scale machine learning. Osdi, volume 16, pages 265-283. Savannah, GA, USA, 2016.
[2] Radhakrishna Achanta, Appu Shaji, Kevin Smith, Aurelien Lucchi, Pascal Fua, and Sabine Süsstrunk. Slic superpixels compared to state-of-the-art superpixel methods. IEEE transactions on pattern analysis and machine intelligence, 34(11):2274-22)) 22822, 2012.
[3] Sushant Agarwal, Shahin Jabbari, Chirag Agarwal, Sohini Upadhyay, Steven Wu, and Himabindu Lakkaraju. Towards the unification and robustness of perturbation and gradi- ent based explanations. In International Conference Machine Learning, pages 110-119. PMLR, 2021.
[4] Chunrong Ai and Edward C Norton. Interaction terms in logit and probit models. Economics letters, 80(1):123 129, 2003.
[5] David Alvarez-Melis and Tommi S. Jaakkola. On the robustness of interpretability methods. CoRR, abs/1800 080490046,.0...... 2018.
[6] Krishnakumar Balasubramanian and Saeed Ghadimi. Zeroth-order nonconvex stochastic op- timization: Handling constraints, high dimensionality, and saddle points. Foundations of Computational Mathematics, pages 1-42, 2022.
[7] Randall Balestriero and Richard Baraniuk. Mad max: Affine spline insights into deep learning, 2018. URL https: // //wrx/aaa.iv..// org/ abs/1805 06576.
[8] Naman Bansal, Chirag Agarwal, and Anh Nguyen. SAM: The Sensitivity of Attribution Methods Hyperparamcters. In Proceedings of the IEEE/CVF Conference on Computer Vision and Pattern Recognition, pages 8673-8683, 2020.
[9] Mahsa Bazzaz and Seth Cooper. Active learning for classifying 2d grid-based level completabil- ity. arXiv preprint arXiv:2309.00 043677 2023.
[10] Andrew Beam and Ben Kompa. Second opinion needed: Communicating uncertainty in medical artificial intelligence. NPJ Digital Medicine, 4, 2021.
[11] Sandeep Bodduluri, Arie Nakhmani, Joseph M. Reinhardt, Carla G. Wilson, Merry-Lynn McDonald, Ramaraju Rudraraju, Byron C. Jaeger, Nirav R. Bhakta, Peter J. Castaldi, Frank C. Sciurba, Chengeui Zhang, Purushotham V. Bangalore, and Surya Bhatt. Deep neural network analyses of spirometry for structural phenotyping of chronic obstructive pulmonary discase. JCI Insight, 5(13), July 2020. doi: 1i17://jiiinnnghh.3378...
[12] Stephen Boyd and Lieven Vandenberghe. Convex optimization. Cambridge university press, 2004.
[13] James Bradbury, Roy Frostig, Peter Hawkins, Matthew James Johnson, Chris Leary, Dougal Maclaurin, George Necula, Adam Paszke, Jake VanderPlas, Skye Wanderman- Milne, and Qiao Zhang. JAX: composable transformations of Python+NumPP programs, 2018. URL http://github.. cww/ google/jax...ooo
[14] Boris Chen, Amir Ziai, Rebecca Tucker, and Yuchen Xie. Match cutting: Finding cuts with smooth visual transitions. 2023 IEEF/CVE Winter Conference Applications of Com- puter Vision (WACV), pages 2114-2124, 2022. URL. httph::tp////c...t //appi www.nticscoolee......... org/
252844493.

---

[15] Noel Codella, Veronica Rotemberg, Philipp Tschandl, M Emre Celebi, Stephen Dusza, David Gutman, Brian Helba, Aadi Kalloo, Konstantinos Liopyris, Michael Marchetti, al. Skin lesion analysis toward melanoma detection 2018: A challenge hosted by the international skin imaging collaboration (isic). arXiv preprint ArXivv9902 0300003368, 2019.
[16] Justin Cosentino, Babak Behsaz, Babak Alipanahi, Zachary R. McCaw, Davin Hill, Tae-Hwi Schwantes-Ain Donghing Lai, Andrew Carroll, Brian D. Hobbs, Michael H. Cho, Cory McLean, and Farhad Hormozdiari. Inference of chronic obstructive pulmonary discase with deep learning on raw spirograms identifies new genetic loci and improves risk models. Nature Genetics, April 2023, ISSN 1546-1718. doi: 10. 1.11008/441588-02--1117272.
[17] Ian Covert, Scott Lundberg, and Su-In Lee. Feature removal a unifying principle for model explanation methods. arXiv preprint 2XXi:2011.006233 2020.
[18] Tianyu Cui, Pekka Marttinen, and Samuel Kaski. Recovering pairwise interactions using neural networks. CoRR, abs/ abss/1901.0861,, 2019. URL httt:://arxiv org/ at////b///// 19901.08661.
[19] Jia Deng, Wei Dong, Richard Socher, Li-Jia Li, Kai Li, and Li Fei-Fei. Imagenet: A large- scale hierarchical image database. In 2009 IEEE conference computer vision and pattern recognition, pages 248-24882 255. leee, 2009.
[20] Ann-Kathrin Dombrowski, Maximilian Alber, Christopher J. Anders, Marcel Ackermann, Klaus-Robert Müller, and Pan Kessel. Explanations can be manipulated and geometry is to blame, 2019. URL https: http:///aaaxvvvorraaaba wwwg/bs////a...eo... hts//bbs/1900 (ttbb///99..0700 07983.
[21] Ann-Kathrin Dombrowski, Christopher J Anders, Klaus-Robert Müller, and Pan Kessel. Towards robust explanations for deep neural networks. Pattern Recognition, 121:108194, 2022.
[22] Charles Dugas, Yoshua Bengio, François Bélisle, Claude Nadeau, and René Garcia. Incor- porating second-order functional knowledge for better option pricing. Advances neural information processing systems, 13, 2000.
[23] Murat A Erdogdu. Newton-stein method: An optimization method for glms via stein's bein's lemma. The Journal of Machine Learning Research, 17(1):7565-7666. 2016.
[24] Gabriel Erion, Joseph D. Janizek, Pascal Sturmfels, Scott M. Lundberg, and Su-In Lee. Im- proving performance of deep learning models with axiomatic attribution priors and expected gradients. Nature Machine Intelligence. 3(7):620-6)) 631, July 2021. ISSN 2522-5839. doi: 10.1038/s42256-02 021-00343-w
[25] Jerome Friedman and Bogdan Popescu. Predictive learning via rule ensembles. The Annals of Applied Statistics, 2, 12 2008. doi: 10.221407-A0AAS448
[26] Bolin Gao and Lacra Pavel. On the properties of the softmax function with application in game theory and reinforcement learning. arXiv preprint arXiv:1704.00805. 2017.
[27] Amirata Ghorbani, Abubakar Abid, and James Zou. Interpretation of Neural Networks Is Fragile. Proceedings of the AAAI Conference on Artificial Intelligence, (100)):661-3688 July 2019. ISSN 2374-3468. doi: 10. 10.1609/aaai, v33101. 333013681.
[28] Xavier Glorot, Antoine Bordes, and Yoshua Bengio. Deep sparse rectifier neural networks. In Proceedings of the fourteenth international conference on artificial intelligence and statistics, pages 315-323. JMLR Workshop and Conference Proceedings, 2011.
[29] Henry Gouk, Fibe Frank, Bernhard Pfahringer, and Michael Cree. Regularisation of neural networks by enforcing lipschitz continuity. Machine Learning, 1 110:393-416, 2021.
[30] Boris Hanin and David Rolnick. Deep relu networks have surprisingly few activation patterns. Advances in neural information processing systems, 32, 2019.
[31] Boris Hanin, Ryan S Jeong, and David Rolnick. Deep relu networks preserve expected length. In International Conference on Learning Representations, 2021.
12

---

[32] Kaiming He. Xiangyu Zhang, Shaoqing Ren, and Jian Sun. Deep residual learning for im- age recognition. In Proceedings of the IEEE Conference on Computer Vision and Pattern Recognition (CVPR), June 2016.
[33] Matthias Hein, Maksym Andriushchenko, and Julian Bitterwolf. Why relu networks yield high-confidence predictions far away from the training data and how to mitigate the problem. 2019.
[34] Dan Hendrycks and Kevin Gimpel. Bridging nonlinearities and stochastic regularizers with gaussian error linear units. 2016.
[35] Davin Hill, Aria Masoomi, Sandesh Ghimire, Max Torop, and Jennifer Dy. Explanation uncertainty with decision boundary awareness. arXiv preprint arXiv:2210.02419 2022.
[36] Davin Hill, Max Torop, Aria Masoomi, Peter Castaldi, Edwin K. Silverman, Sandeep Bodduluri, Surya P Bhatt, Taedong Yun, Cory Y McLean, Farhad Hormozdiari, Jennifer Dy. Michael Cho, and Brian D Hobbs. Deep learning utilizing suboptimal spirometry data to improve lung function and mortality prediction in the uk biobank. medRxiv, 2023, doi: 10.1 10.1101/20 04, 28.23289178. URL https: //ww: medrxiv. org/ //t////n/eaarrt////////////3//// 28. 23289178.
.4.
[37] Cheng-Yi Yu Hsieh, Chih-Kuan Yeh, Xuanqing Liu, Pradeep Ravikumar, Seungyeon Kim, Sanjiv Kumar, and Cho-Jui Hsieh. Evaluations and Methods for Explanation through Robustness Analysis. arXiv:2006.00442 [cs, stat], April 2021.
[38] Zihan Huang, Charles Low, Menggiu Teng, Hongyi Zhang, Daniel E. Ho, Mark S. Krass, and Matthias Grabmair. Context-aware legal citation recommendation using deep learning. CoRR, abs/2106.10776, 2021. URL https: / org/ hbtp////////....// 10776.
[39] Joseph D Janizek, Pascal Sturmfels, and Su-In Lee. Explaining explanations: Axiomatic feature interactions for deep networks. J. Mach. Learn. Res., 22:104-1: 2021.
[40] David P. Johns, Julia A. E. Walters, and E. Haydn Walters. Diagmosis and early detection of COPD using spirometry. Journal of thoracic disease, 6(11)1557,, November 2014. ISSN 2072-1439 2077-6624. doi: (..3978/..ssn.22722.. 143.02--144...11.08..1
[41] Zulqarnain Khan, Aria Masoomi, Davin Hill, and Jennifer Dy. Analyzing the effects of classifier lipschitzness on explainers. arXiv preprint (.Kiv:2206.1248.1 2022.
[42] William Knauth. The self-simplifying machine: Exploiting the structure of piecewise linear neural networks to create interpretable models. arXiv preprint arXiv:2012.01299 2020.
[43] Alex Krizhevsky, Geoffrey Hinton, et al. Learning multiple layers of features from tiny images. 2009.
[44] Jeffrey Larson, Matt Menickelly, and Stefan M. Wild. Derivative-fres optimization methods. Acta Numerica, 28:287-400 2019. doi: 10.1017/509624929190000 19000000.
[45] Yann LeCun and Corinna Cortes. MNIST handwritten digit database. 2010. URL http:
/ //yann. lecun. com/ hamm//xxb//mnis//////
[46] G. Leoni. A First Course in Sobolev Spaces. Graduate studies in mathematics. American Mathematical Soc., 2009. ISBN 978-0-8218-8415-7 URL https://books.. google com/ books??id=W3RL.WmYYOR
[47] Samuel Lerman, Charles Venuto, Henry Kautz, and Chenliang Xu. Explaining local, global, and higher-order interactions in deep learning. In Proceedings of the IEEE/CVF International Conference on Computer Vision, pages 1224-1233, 2021.
[48] Cong Han Lim, Raquel Urtasun, and Ersin Yumer. Hierarchical verification for adversarial robustness. Hal Daumé III and Aarti Singh, editors, Proceedings of the 37th International Conference Machine Learning, volume 119 of Proceedings of Machine Learning Research, pages 6072-6082, PMLR, 13-18 Jul 2020. URL https: hwtph:tttpppcceeeeee h/proeeed//cccocc.g... lim20b. html.
.aweeeddiinnm..........mm.m.mere mlr. prrss//ve1/
13

---

[49] Wu Lin, Mohammad Emtiyaz Khan, and Mark Schmidt. Stein's lemma for the reparameteri- zation trick with exponential family mixtures, 2019. URL htttpp:/aarxiv.orggg.// abgg//////11 13398.
[50] Zachary C Lipton. The mythos of model interpretability: In machine learning, the concept of interpretability is both important and slippery. Queue, 16(3):3 2018.
[51] Jun S Liu. Siegel's formula via stein's identities. Statistics & probability letters, 21(3), 1994-10. ISSN 0167-7152.
[52] Ilya Loshchilov and Frank Hutter. Sgdr: Stochastic gradient descent with warm restarts. arXiv preprint arXiv: 03983, 2016.
[53] Scott M Lundberg and Su-In Lee. A unified approach to interpreting model predictions. Advances in neural information processing systems, 30, 2017.
[54] Scott Lundberg, Gabriel Erion, Hugh Chen, Alex DeGrave, Jordan Prutkin, Bala Nair, Ronit Katz, Jonathan Himmelfarb, Nisha Bansal, and Su-In Lee. From local explanations to global understanding with explainable ai for trees. Nature machine intelligence, 2(1):56-67, 2020.
[55] Andrew Z Luo, Eric Whitmire, James W Stout, Drew Martenson, and Shwetak Patel. Automatic characterization of user errors in spirometry. In 2017 39th Annual International Conference of the IEEE Engineering in Medicine and Biology Society (EMBC), pages 4239-4242. IEEE, 2017.
[56] Aria Masoomi, Chich Wu, Tingting Zhao, Zifeng Wang, Peter Castaldi, and Jennifer Dy. Instance-win wise feature grouping. Advances in Neural Information Processing Systems, 33: 13374 13386, 2020.
[57] Aria Masoomi, Davin Hill, Zhonghui Xu, Craig Hersh, Edwin Silverman, Peter J Castaldi, Stratis loannidis, and Jennifer Dy. Explanations of black-box models based on directional feature interactions. International Conference on Learning Representations, 2021.
[58] Guido F Montufar, Razvan Pascanu, Kyunghyun Cho, and Yoshua Bengio. On the number of linear regions of deep neural networks. Advances in neural information processing systems, 27, 2014.
[59] V.C. Moore. Spirometry: step by step. Breathe, 8(3):232-240) 2012. ISSN 1810-6838. doi: 10. 1183/20734735.0022 0021111.. URL https: ers jrsjonralss com/content///3/222.
[60] Yurüi Nesterov and Vladimir Spokoiny, Random gradient-free minimization of convex functions. Foundations of Computational Mathematics, 17:527-566, 2017.
[61] Daniel Omeiza, Skyler Speakman, Celia Cintas, and Komminist Weldermariam. Smooth grad- cam++: An enhanced inference level visualization technique for deep convolutional neural network models. arXiv preprint 1908 0/244, 2019.
[62] Adam Paszke, Sam Gross, Soumith Chintala, Gregory Chanan, Edward Yang, Zachary De Zeming Lin, Alban Desmaison, Luca Antiga, and Adam Lerer. Automatic diffcrentiation in pytorch. 2017.
[63] Hanchuan Peng, Fuhui Long, and Chris Ding. Feature selection based on mutual information: Criteria of max- max-dependency, max-relevance, and min-redundancy. IEEE Trans. Pattern Anal. Mach. Intell., 27(8): 1238, aug 2005. ISSN 0162-8828, doi: 10. 20.1009/AM....00.115555 URL https / /doi. org/10 10. 1..0..1.0 1119/ 2005. 159.
[64] Maithra Raghu, Ben Poole, Jon Kleinberg, Surya Ganguli, and Jascha Sohl-Dickstein. On the expressive power of deep neural networks. In international conference on machine learning, pages 2847-2854. PMLR, 2017.
[65] Prajit Ramachandran, Barret Zoph, and Quoe V Le. Searching for activation functions, arXiv preprint arXiv:1710.05941, 2017.
14

---

[66] Marco Tulio Ribeiro, Sameer Singh, and Carlos Guestrin. "why should trust you?" explaining the predictions of any classifier. In Proceedings of the 22nd ACM SIGKDD international conference on knowledge discovery and data mining, pages 1135-1144, 2016.
[67] Abhijit Guha Roy, Jie Ren, Shekoofeh Azizi, Aaron Loh, Vivek Natarajan, Basil Mustafa, Nick Pawlowski, Jan Freyberg. Yuan Liu, Zach Beaver, et al. Does your dermatology classifier know what it doesn't tt know? detecting the long-tail of unseen conditions. Medical Image Analysis, 75: 102274, 2022.
[68] H. Royden and P. Fitzpatrick. Real Analysis. Pearson Modern Classics for Advanced Mathe- matics Series. Pearson, 2017. ISBN 34689944 URL https://books: google. com/ books?7dd=7GP3MAAAAA.
[69] Patrick Schwab and Walter Karlen. CXPlain: Causal explanations for model interpretation under uncertainty. In H. Wallach, H. Larochelle, A. Beygelzimer, F. dAlché-Buc, E. Fox, and R. Garnett, editors, Advances in Neural Information Processing Systems, volume 32. Curran Associates, Inc., 2019.
[70] Ramprasaath R. Selvaraju, Abhishek Das, Ramakrishna Vedantam, Michael Cogswell, Devi Parikh, and Dhruv Batra. Grad-cam: Why did you say that? visual explanations from deep networks via gradient-based localization. CoRR, abs/1610.02391, 2016. URL http: //axiv.. org/abs /1100. 023991.
[71] Thiago Serra, Christian Tjandraatmadja, and Srikumar Ramalingam. Bounding and counting linear regions of deep neural networks. In International Conference on Machine Learning. pages 4558 4566. PMLR, 2018.
[72] Karen Simonyan and Andrew Zisserman. Very deep convolutional networks for large-scale image recognition. arXiv preprint 1409v 2014.
[73] Karen Simonyan, Andrea Vedaldi, and Andrew Zisserman. Deep inside convolutional networks: Visualising image classification models and saliency maps. CoRR, 1bl.644,,, 2014.
[74] Sahil Singla, Eric Wallace, Shi Feng, and Soheil Feizi Understanding impacts of high-order loss approximations and features in deep learning interpretation. In International Conference on Machine Learning, pages 5848-5856. PMLR, 2019.
[75] Dylan Slack, Sophie Hilgard, Sameer Singh, and Himabindu Lakkaraju. Reliable Post hoc Explanations Modeling Uncertainty in Explainability. In Neural Information Processing Systems (NeurIPS), 2021.
[76] Daniel Smilkov, Nikhil Thorat, Been Kim, Fernanda Viégas, and Martin Wattenberg. Smooth- grad: removing noise by adding noise. arXiv preprint 13Xxv:1706.03888 2017.
[77] Daria Sorokina, Rich Caruana, Mirek Riedewald, and Daniel Fink. Detecting statistical inter- actions with additive groves of trees. In Proceedings of the 25th international conference on Machine learning. pages 1000-1007, 2008
[78] Charles M Stein. Estimation of the mean of a multivariate normal distribution. The annals of Statistics, pages 1 1135-1151, 1981.
[79] Cathie Sudlow, John Gallacher, Naomi Allen, Valerie Beral, Paul Burton, John Danesh, Paul Downey, Paul Elliott, Jane Green, Martin Landray, Bette Liu, Paul Matthews, Giok Ong, Jill Pell, Alan Silman, Alan Young. Tim Sprosen, Tim Peakman, and Rory Collins. Uk biobank: An open access resource for identifying the causes of a wide range of complex diseases of middle and old age. PLOS Medicine, 12(3):1-10, 03 2015. doi: 10. 10.13711 1/jor//l pmed. mimedd.101171 779. URL https: http:///oooooo/...////
/ pmed. 1001779.
[80] Adam Summerville, Sam Snodgrass, Matthew Guzdial, Christoffer Holmgård, Amy K Hoover, Aaron Isaksen, Andy Nealen, and Julian Togelius. Procedural content generation via machine learning (pegml). IEEE Transactions on Games, 10(3):557-770, 2018
[81] Mukund Sundararajan, Ankur Taly, and Qiqi Yan. Axiomatic attribution for deep networks. In International conference machine learning, pages 3319-3328. PMLR, 2017.
15

---

[82] Mukund Sundararajan, Kedar Dhamdhere, and Ashish Agarwal. The Shapley Taylor Interaction Index. In Proceedings of the 37th International Conference on Machine Learning, pages 9268. PMLR, November 2020.
[83] Tijmen Tieleman, Geoffrey Hinton, et al. Lecture 6.5-rmsprop: Divide the gradient by running average of its recent magnitude. COURSERA: Neural networks for machine learning, 4(2):
26-31. 2012.
[84] Michael Tsang, Dehua Cheng, and Yan Liu. Detecting statistical interactions from neural network weights. In International Conference on Learning Representations, 2018.
[85] Michael Tsang, Dehua Cheng, Hanpeng Liu, Xue Feng, Eric Zhou, and Yan Liu. Feature interaction interpretability: A case for explaining odfrecommennation systems via neural interaction detection. In International Conference on Learning Representations, 2019.
[86] Michael Tsang, Sirisha Rambhatla, and Yan Liu. How does this interaction affect me? inter- pretable attribution for feature interactions. Advances neural information processing systems, 33:6147-6159, 2020.
[87] Roman Vershynin. Introduction to the non-asymptotic analysis of random matrices. arXiv preprint 101.:71 2010.
[88] Roman Vershynin. #ligh-dimensional probability: An introduction with applications in data science, volume 47. Cambridge university press, 2018
[89] Zifan Wang, Haofan Wang, Shakul Ramkumar, Piotr Mardziel, Matt Fredrikson, and Anupam Datta. Smoothed geometry for robust attribution. Advances in Neural Information Processing Systems, 33:13623-11 13634,, 2020.
[90] Han Xiao, Kashif Rasul, and Roland Vollgraf. Fashion-mnist: a novel image dataset for benchmarking machine learning algorithms. CoRR, abs/1708.07747, 2017. URL http://
wexiv..rg/abee
07747..
[91] Taedong Yun, Justin Cosentino, Babak Behsaz, Zachary McCaw, Davin Hill, Robert Luben, Dongbing Lai, John Bates, Howard Yang, Tae-Hwi Schwantes-An. Anthony P Khawaja, An- drew Carroll, Brian D Hobbs, Michael Cho, Cory Y McLean, and Farhad Hormozdiari, Unsupervised representation learning improves genomic discovery for lung function and respiratory disease prediction. medRxiv, 2023. doi: 10. 11111//222..4...2233333333 https://vet http:://///mmmmdziiv
URL 04 28. 23289285.
[92] Hao Zhang, Yichen Xie, Longjie Zheng, Die Zhang. and Quanshi Zhang. Interpreting multi- variate shapley interactions in dnns. In The AAAI Conference on Artificial Intelligence (AAAI), 2021.
[93] Bolei Zhou, Aditya Khosla, Agata Lapedriza, Aude Oliva, and Antonio Torralba. Learning deep features for discriminative localization. In Proceedings of the IEEE conference computer vision and pattern recognition, pages 2921-2929, 2016.
[94] Jingyi Zhu. Hessian estimation via stein's identity in black-box problems. arXiv preprint arXiv:2104.01317, 2021.
16

---

A Societal Impacts
Machine learning models have become increasingly pervasive in society, from medicine [10, 67, 15] and law [38] to entertainmen([4, 14, 9, 803. Therefore, it is important that users of the technology understand the factors underlying model predictions. To this end we propose SmoothHess for quantifying the feature interactions affecting model output. Potential applications of our method are widespread; SmoothHess may be used to find interactions influencing model to predict whether customer will default a credit loan or if patient has melanoma. The deeper understanding of the model gleaned from SmoothHess may used to improve decision making. For instance, doctor may notice that the SmootHHess feature interactions between the pixels in image of lesion don't "make sense", indicating that the model should not be trusted and that a more granular human assessment is required. However, such applications are highly sensitive and the cost of inaccurate predictions, or, in this case, misinterpreted attributions, can be high. An inaccurate interpretation can instill unwarranted trust, mistrust, in user. In the most extreme cases this may lead to sub-optimal decision making such as the confident denial of credit loan to a trustworthy customer misdiagnosis of benign lesion malignant.
In addition, it is important to note the challenges related to ensuring robust, stable, and trustworthy explanations. In particular, recent works have uncovered issues related to the sensitivity of the explainer to small changes the input [5, 41], adversarial attacks [27, 20, 37], or hyperparameter tuning [8]. Methods have been proposed that attempt to quantify explanation uncertainty [75, 69, 35], however further challenges remain. Thus, as with all methods for explaining machine learning model predictions, we recommend that SmoothHess is used in tandem with the careful consideration of domain experts, who are best equipped to interpret interactions the context of their field.
B Proof of Proposition 1
B.1 Preliminary
We make use of a Lemma from Lin et al. [49] in our proof below, which we relate here:
Lemma 3. (Lin et al. [49]) Denote xo Ed, Docally-Lipschite continuous function h(2) R2d R.
covariance matrix E RRdxd random vector IRd distributed from NZ E). Ef [I(()]]] < 20, then
E ((z ((z-x))(> z0) = xo) (2- V.)a)""
(13)
B.2 Main Result
Proposition 1. Given xo E R24, L-Lipschitz continuous functiinn g Rd R. covariance matrix Σ E Raxd and random vector 5 RR distributed from 5 [0, with density function gr. Rdd R. then
EA E V. + 6)) - (x))) Vohg,((0o).
where * denotes convolution.
Proof. We define random vector = = Fn + 8 E Rd distributed from Nixx: E), and for which
we denote the density function py. We begin by showing [[]]]] Esslg(z + 6)]l < co. We are given that is L-Lipschitz, i.e. Va, d E Idd
Ig(a) g(b) S Llla b112
(14)
for some L >>0 0.
Now, for any fixed 8 E ddd we have
ineq.
Eq.14
I0( + 5) 1g(=0) + (g(xo 8) - 9(x))) - I9(x0)] + Iø( + 6) g(zo)l
(15a)
1g(0))) + -| + =oll2 = (g(2)) ) + L||6||2-
(15b)
Thus, we may bound the expectation Es ( + 6]]]:
Esllg000 + +]]] < [s(lo(=))) + L [| = g(æa) + [Es[11112] < cc.
(16)
17

---

Here ((20) < oc as it is a constant, and it may be seen that LEg [|]]||]] < co using simple change of variables. Defining 3 Σ' d ~ 10, one may write Es [ | As
straightforward consequence of Cauchy-Schwarz for any fixed 3, one may write
د> B11 - ||> |||||||
(17)
where ||-|| denotes the Frohenius norm. Noting that || Xa Eq. (17) to see
I [(dd+((())))) 1)/2)
Eg 1 A|la] s [> || | مل و|
2
< 000 (18)
T(d/2)
where [() denotes the Gamma function.
Next, we move the Hessian operator inside the integral.
V2 + gx)) = v2 /
g(2)g=(z zo)dz /
2ZERe g(z) [(= - -z0)d==
(19a)
ZERS
/ V mn)dz 1
(V xx)((= -5 ---dd==
(19b)
EERR
JEERU
f ((=)((V<( =00 ))2 (3 Io) + (V~qz(>~ - zo))(3 20) "00 "d)
(19c)
ZERU
ZERR ((=)(-12 2 mo) + (V xo))(z 20)² 2 dd==
(19d)
f g(2)(- 12 mm) + 93(2 20)=-11 (2 - (2 - I
(19e)
"ERR
[ g(2)91(2 monr + Σー ( (z x0)(2 =00²² " ))z==
(19f)
JEER
[ (1) -- (2
(2 m
(px(2) qr(2 xn))
(19g)
SSERd
B.]g(=) 13 + -Σー - mn) xo) T =
(19h)
E (-I (2 1 x0) - -00) -
(191)
(-E+
xo) (1))))) 1)
(19j)
Ex ((2 (====o)(= z0)² -ΣソΣ [(]]]
(19k)
Lemma 3 may be applied as g is Lipschitz, and thus locally-Lipschizz, and E=(lo(=)) < co, yielding
E (= =00²² -EEE- E2 (z0))_<
(20)
Using change of variables from = to xi + 8 we write
+ 6 which, when combined with Eq. (19) and Eq. (00), completes the proof
E E2 - x0) V [=(())]]]] =
(21)
22 * (zo) Es E +))
(22)
C Proof of Theorem 1
We begin by establishing a result expressing the eigenvalues of the symmetrization of a rank matrix in closed form:
Lemma 4. Given F, E Rd denote A = =xy" + yzu E Rddd, The following facts hold:
1. Matrix A will have the following eigenvalues
It has d - 2 eeigenvalues equal 0
18

---

The other two eigenvalues are denoted by x (A) and X (A) will have the following form:
X (A) - x"y + |
(23)
2. (A) and X (A) are non-negative and non-positive respectively.
3. Given and y are sampled from sub-gaussian distributions, then A (A) and X (A) are sub-exponential random variables.
Proof. We divide the proof into three sections
1. We have (A)) 2, A E So- Therefore EQ [e.lee].../eel Rddd s.t. Q²Q
QQ" = T, the column vectors e; E RR are orthonormal and
A= QAQ" Q"AQ=
(24a)
(24b)
where A = diag([A1, A2, ... 0)) E RRdxd As eigenvalues are invariant to change of basis, has eigenvalues A₁ and X2 and the other 2 eigenvalues are equal to O.
It can be seen that span( {.x, y = span(fe1 e2 As span((e), e2 = C(A) this can be
shown by proving span {z, y}) = C(A). We know 3zg E Rd : 2x 1 x and Bzy E Id
By 1 y. Thus we have
Az = (any T + yx "zz")zz-> = xy 2x + J=
2 x(y xx)
(25a)
Azy (ay² + yz²))ay )zy my By + yx By y(x zny)
(25b)
From the above, have x * and thus span(fr: y/) = C(A). If rank A) = we have x * 0, + 0 and still span( (z. v/) C(A). If (A))(A) 0 we have x y = 0 and clearly span (fz, y/) {0} C(A).
that 2;, y E C(A). We know rank(A) S 2. If rank(AA 2
As {e are orthonormal, we have that e, y" ej 0, Vj 2. We define
ã = Q² x ae e (T) = (eT o x, 0, ................................ 0) E Rd.
(26a)
ÿ Q" (ef o v, Tv. U ...
. 0) Rd
(26b)
We define Q2 = [e1 E ddx2 and
= = x (ef e E TR?
(27a)
ÿ ÿ= T y = (eT (eTu.e5 v ER²
(27b)
We first show that || and l are equal to ||0||| and l|llle respectively:
||ا|| = l1Q 21|2 = ((QTx) x) = V=TQQT==
|||1|
(28a)
||7||2 ||Q" ylla = V ()TT) V VTQQQT = V y ||yll2
(28b)
Next, we use Eq. (28) to show that || and |lūll2 are equal to ||ا||| and llllla respectively:
|لا| = V (eTx))+ = | = ||0|||
(29a)
|b - (cFu) + (cTv)² 10 = || la
(29b)
We show an equality between inner products 3 = x" y:
T ア = n" Ty = ="y
+ +e
= X" " (Q" x) Q" Q²= = -QQ"==
(30a)
(30b)
19

---

Now, one may write,
AQ = diag( (A1, A2)) 0 (xy" + yx ]Q2 yr²))2 = diag( (1, >z))
(31a)
(31b)
2y Q2 + Q2 12yyz" Q2 diag , A2))
(31c)
an" gā" diag (0 X21))
(31d)
The following facts hold as result Eq. (31d):
\ 2.
A2 2T:222
T182 + T28l1 0
Now, show
A = ," | Which can be derived as such:
A2 3 - |
(32)
q"ツ + | = Z101
T2ÿ2 +
55
(33)
= ē1 T202 +
(Z1ğ1 و2 2
3 01 + おхうう + ت1
A A2
where the second equality comes from the following:
(Σ T2222)² (ā1
27 22, (ū2う2)²
(Z101)² 27, + (i2222) + (Z192 22011²
(34)
( + (1 + 52
Proving that Eq. (32) holds. Finally, we combine Eq. (29) and Eq. (30) with Eq. (32) to
express the cigenvalues of A as:
== x=y y + ||,ا||/|ا||
2 = xJ - |0||2|111z
(35)
which we use denote X (A) A ((A) = A2.
2. As Jx yl < | it follows that x+)) xTy | 2 0 and X (A)
x y | llvll $ 0.
3. We denote D = (1....,d d} can be seen that 2.3// is sub-cxponential Vi e as sub-gaussian times a uub-gaussian is sub exponential. Thus, it follows that
Tili is a sub-exponential random variable,
(36)
that xf and y are sub-exponential as the square of a sub-gausee sub-gaussian is sub-exponential. As the
the sum of sub-exponential random variables is sub-exponential. Further
sum of sub-exponentials is sub-exponential we have that E 2f; Ed vē are both sub- exponential random variables. As the square root of a sub-exponential is sub-gaussian we have that
||0|||
x is sub-gaussian random variable
(37a)
|lvl|2 vf is sub-gaussian random variable
(37b)
As a sub-gaussian times sub-gaussian is sub-exponential, from Eq. 37 we have that ||ا||||||||| is sub-exponential random variable.
(38)
Now we see from Eq. and 36 Eq. 38 that both x+( and X >-(A) are the sum of sub- exponential random variables and thus are sub-exponential.
20

---

We now use the result of Lemma 4 to prove the sample complexity bounds for SmootHHess in Theorem 1:
Theorem 1. Let f
R R be a piece-wise linear function over finite partition of R.. Let Do E Rd, and denote {diJiz1. set of n ii.d random vectors in R distributed
from o; ~ E). Given H. as in Eq. (8), for any fixed E,y E (0, 1), given then
2 A [maxx Vd+ ct log 2), (C vd+ v
_ log 2]
[| H1|2 >8) S/
where H VVIff (v) C+, 20,, > covariance E and qE : R + R the density function of (0)
are constants depending on the function f and
Proof. As To, and E are fixed, refer to Hn Io, 2) Hn for brevity. We denote D (l....,dl d3. We begin by explicitly expressing our estimator H, in terms of d; and V (2o 6:).
From Eq. (8) we have
a H
E
8 V=, + 8.)"'
(39a)
n
1=1
Hn
H
H
(39b)
2
1 E & V [ + 8,)
11
+ 2 1=1
Σ
f + 6,)o,' x-1)
(39c)
2 n
1=1
Now, we show the convergence of our estimator An: Lemma 5, limm-xx Hn H
Proof. From Proposition it is clear to see that limm H HI:
lim H lim
Sr-'o,IV2 + +6)) EE(E- of -
+ +)17 Prop
H. (40)
n-000 n-000
1=1
Next, we show it is also the case that lim, 11 H::
lim H. lim
V
(zo +8.)64--1-
(41a)
lim 6 V + 5.) = lim
6, V f(a + 60)
(41b)
n1000
1=1
1L + 1)) Prop HT H
(Symmetry of Hessian)
(41c)
Now, as we have lim H liim HoT
H, it follows that
lim H. lim 1 H lim
1 HoT
1
H
--
H
(42)
2
200
2
2
We establish the following notation be used below: given fixed vector 8 Id one may construct matrix As kd by:
As 1 (E- + 6)1²) + 1 (V) 8-8- Rdxd
(43)
21

---

It can be seen from Eq. 39 and Lemma 5 that H and Hn, may be expressed in terms of matrices As and Ag.:
H = R.s [A6],
Hn
A6.
(44)
Next, we establish that the random vectors Σ- and (2o + 6) are sub- sub-gaussiam: Lemma 6. The random vectors 1 and V fxo + 8) are sub sub-gaussian.
Proof. As E-15 is Gaussian it is sub-gaussian. Now, we show that V f100 8) is sub-gaussian random-vector We have that is piecewise-linear over a partition of IR with finite cardinality L. Let us denote this partition as Q = (/ Qi C IR", where, when restricted to given Q E Q we have
f - V + Aq
(45)
where Vi22 R2, Aq E R are the affine coefficients associated with the region Q. Then it is the case that V Rd R is bounded function, where, aside from set of measure 0, M Rd (the boundaries of regions Q) where Vf is not defined, one has
|V)(2=) s max Vollz. Var E M.
(46)
QEQ
Thus, Vf(o + 5) is a bounded random vector and therefore is sub-gaussian.
Given the operators X A Raxd + R as defined in the statement of Lemma 4 and fixed vector Rd, we denote 2 == X+ (As), := X (As) and the corresponding unit eigenvectors as vgs E Rd and V/3 E IRd respectively. We denote random vectors w w1 E IRd by
w
Wg
A U/8
(47)
where v A and of are random variable random vector pair coming from the same 8. An immediate consequence of Lemma 6 is that w and wig are sub-gaussian random vectors:
Lemma 7. w and wg are sub-gaussian random-vectors
Proof. Using Lemma 4(3), Lemma 4(2) and Lemma we see that that A and -Xī are non-negatin nen-negative sub-exponential random variables and thus that / /x and V-x are sub-gaussian random variables. We may say w is sub-gaussiaa random vector if (wT :) is sub-gaussian random variable Vz e Rdd Let us fix arbitrary z E Rdd As A is sub-gaussian we have
3K > 0st. P( A 2 f) 2exp(tt/1)))
(48a)
Now, Vt 2 0
(0 ; =) 2 P(ly 13 (0f =) 2)
(49a)
C-S. || 1l
F(1 At 2 t/ 2)1)
P(
1 t/1+111222 22exp(-(P/1l)ff)
(49b)
Thus defining K(a) = K2 ||21|2 we see that 2 0
P(Kkw 2)|2 2 1) < exp((t2/(K)))
(50)
Thus (wz 2) sub-gaussian for arbitrary E Rd. Therefore, w a sub-gaussian random-vector. The same argument holds to show that wiz is a sub-gaussian random vector.
22

---

Given any fixed 6, it can be seen that
As = w/ 4T
w W8
(51)
Thus one may wrrit H and H, from Eq. 44 as:
1
n
1Es [As] (w3 w/" w wg ". Hn
A6:
w t +T b. 6
-, wb
(52)
11
Because Vi=1 are i.i.d. random vectors and w " are fully determined by 6,, it follows that {wg, Yi= and {w5 /i are both sets of i..d. random vectors.
We denote the following:
H Es(ug w '].
H Eg [wg W/ 7].
(53a)
1
1
H1
W w +T
WK, W6: Hr
(53b)
i=1
can be seen from the RHS of Eq. (52) that
Hn = 1 11,
(54)
We now aim decompose H in terms of 11+ 11,,"""" in order to derive separate concentration bounds. To this end, we prove the following lemma:
Lemma 8. H and H exist.
Proof. Let us consider the random matrix wz "0 E Rdxd, As (1x5) )x E TR is sub-gaussian Vk D,
the element (wg " )yg E R is sub-exponential as the product of two sub-gaussian's Vi,j E D. Thus, [s[(w] w +T) )ay] exists Vi,j E The same argument can be made to show Eg (wn wg )ol
exists Vi D.
In light of Lemma 8, the LHS of Eq. (52) may be decomposed as
= Eg (wo w wj w 7 = Es [ w 7 - Es [wy w 7 H
(55)
Before deriving separate concentration bounds on Ht and fr we note that
lim म Es(w] 10,"" <-0
lim म. Eslu wy -1. 200
(56)
Finally, we bound the deviation of Āt and i1 from their expectations. Let us fix 8, E (0, 1]. From Eq. (55), and the fact that {w ) is a set of i.i.d. sub-gaussian random vectors, Theorem 3.39, Remark 3.40 Vershynin [87] may be applied, yielding: Vt 2
P
(
H H+1l2 max(s;t (Et))) < 22exp(-ctt)
(57)
where et = - yd Vx
and C+ c4 constants depending the sub-gaussian norm of Plugging into Eq. (57), we get
Vn
w Let select t
lor(4/>) ol
log(4/7)
log(1/7) ct
Vd
Vd
P | - 11 > +
ct
(0+
)2 -
(58)
Vn
Vn
n
23

---

Let us consider n+
4 (c* Vd + V /og(4//132
(c+
One may see that
v /1og(1/1)
+ 2t = Vd
at
(59a)
Vn+
Vitt
VV log(4/7)
EVd
C
(59b)
2(C Va +
(18(4/11)) at
2(C VC+ / (18(44/7)
C+Vd+ + log(4/7)
(59c)
C+vd+
log(4/7)
2
at
Thus, given n n+ one has ot = max($ (5))) = 3. because & << 1, and that
H | -
E ? 2
(60)
2
In fact, because ot is monotonically decreasing in given 2 n+ Eq. (60) holds.
The same logic above may be used to show that there exists constants C c depending on the sub- sab-gaussian norm of w'y such that, given n 2 n = (C vd log(4/m))2 one has
8 2 Ar H 1l2
?
(61)
Finally, combine the two bounds from Eq. (61) and Eq. (60). Given Z mC* Vd
log(4/7) ) - (C- + lmg(4//) )2) we have
F(I H > 2) P(] 11 - 11 + 11 =))=))=
(Eq. (55) (54)) (62a)
P(I H+) (H, H H1 =)
(62b)
F(I H*) + 11t, - H )
(A Ineq.) (62c)
F(I H+ 122 >
8
U || - H (
(62d)
2
=
5
- - 2)
-
(Union bound) (62e)
2
2 +
(Eq, (60)Eq. (61)) (62f)
2
D
Implementation Details
D.1 Quadratic Optimization
Given a function f Rdd - IR, point no E Rd gradient Hessian pair G E Rd H E Rdxd and a
magnitude constraint > O. we aim to solve the following optimization:
min AERd
f(zo) + GT - mp) + 2 (A (A- 200) - 2o), s.t.
I|A - xolla Ss, (63)
as (zo) constant, the problem above equivalent
min SER4
G H6, s.t. ||01|2 I E,
(64)
where we have replaced To. which can be interpreted as the output after the attack with 6 = A Zo. the attack vector itself.
24

---

The optimization problem in Eq. (64) is non convex H is not guaranteed to be positive semi- definite. However, as Slater's constraint qualification is satisfied, i.c. 26 E IR s.t. ||ö]|2 < E, Eq. (64) may be solved exactly(12]. Specifically, the solution may be obtained by solving an equivalent convex optimization problem:
1
min tr( HX) G VER",XES" 2
(65a)
s.l. (r(x) 8² 0,
(65b)
|x,m7 To =1
(65c)
where 5d denotes the the set of symmetric matrices and ł indicates the block matrix 1x, 1] (dd+1)x(d+1) is constrained to be positive semi-definite
However, the optimization in Eq. (65) is expensive solve when >>>>p 0 as there are ( variables. For instance, MNIST and FMNIST have dl x 6.0 10" and CIFARIO R10 has d² se 10 Thus, before converting Eq. (64) into Eq. (65), elect reduce the dimension of the optimization problem k E K << d.
Let us consider the cigendecomposition QAQT= = H. Here the columns of Q Rd are orthonor- mal eigenvectors of H. Given the d cigenvalues {A
sorted such that i < = | 2 Ayl. we have A diag(A1,-...) . _d)). We remove the last d k columns from Q and d - k columns and rows from to construct Q E R k A R Thus, we have low-rank approximation of H:
H = QAQ".
(66)
Thus, Eq. (64) is approximately equivalent to another optimization which uses this low-rank approxi- mation for H:
G 1 o QAQO' s.t. ||0122 Se.
min SERd
(67)
Defining 5 = TT E RA Eq. (67) is approximately equivalent to
G²O² 2 A8, stt. ||ā]|e S,
min
(68)
SER*
where the constraint is simplified to ||ō||2 from | S c
|l
= (5 Q" Q5) (8"Q" (8' = =118|2-
(69)
Finally, 8* E R the optimal solution Eq. (68), is projected back to Rd yielding an approximate solution to Eq. (67):
8 Q Rd.
(70)
Choosing k: The choice of k E N is determined using threshold hyperparameter 7 E (0, 1]. Given T. is chosen to be the smallest number of (sorted) sigenvalneseeigenvettoss which account for proportion of the total eigenvalue magnitude that is at least T:
K' E
k = argmin
(k' E N:
(71)
1=1
i=1
For the values of T used for each dataset see App. E
D.1.1 Similarities and Differences with CASO/CAFO
Our quadratic optimization closely connected the CAFO/CASO explanation vectors proposed by Singla et al. [74]. Both methods use optimizations to minimize value outputted by the network, and the CASO method also uses second-order Taylor expansion in their objective. Additionally, the proposed Smooth-CASO method is equivalent to attack using SmoothGrad and SmoothHess on the cross-entropy loss function (when their regularizers are set to 0), which admits higher order derivatives.
25

---

However, there are some key differences. First, our optimization is meant for use with arbitrary network outputs as opposed to just the loss, as the case with CAFO/CASO. While the vanilla loss Hessian modeled by Singla et al. [74] is proven to be positive semi-definite, this is not the case for arbitrary functions which one may wish to attack. Thus methods to optimize non-convex quadratic objectives, such as those outlined in App. 1, are required. Second, our goal is not to generate an explanatory feature importance vector for analysis, as an important motivation for CAFO/CASO, but to assess the quality of the gradient Hessian pair used to attack the function. For this reason, we do not use sparsity constraints such as in Singla et al. [74], which are part meant to improve the interpretability of CAFOICASO as explainers. Last stress that the techniques used to find the cross-entropy loss Hessian for CASO and Smooth-CASO cannot be used for internal neurons, logits or regression valued output in ReLU networks, due to their seecewsse-linarrin
D.2 SmoothHess
The SmoothHess estimation procedure, and the amortization with SmoothGrad estimation, is pre- sented below in Algorithm I. While empirically SmoothHess estimator converges, have additionally found that reflecting each point 6, in the perturbation set {6. { about the origin to create an augmented perturbation set {6. U {- before estimation can result in faster per-sample convergence.
Algorithm Joint SmoothHess and SroothGrad Estimation
Input: Sample of interest r E R Neural network indexed to output scalar of interest f R - R, Covariance Σ E Raxd Batch size for gradient oracle calls nı E N, Number of batches n2 E N.
Output an estimate SmoothHess, G an estimate SmoothGrad
A, c torch.zerosdd.d), torch.zeros(d) E-1 + torch.inverse(E))
\\ H E 1x G E Ied
|1 1 E Raxd If 5. diagonal (dd) Else (dd)
110 (m + d)), O(W) is complexity of one forward pass through f for
122 do
8 + torch. 0, Σ)
115 Rixd
vf(= 00 + torch f 6)
Vf e O(n,
Å. torch.matmul(6, vf(x 6).T)
d E Raxd O(nid')
i < f + A. m112
\1 (n)dd)
â + G Vf(z ( + summ 0) n1722
((ndd)
end
A torch torch..atmul(2-1_00 H)
\\ If diagonal O(d?), Else (dd)
A - A H.T
Retum A, G
D.3 Alternative Covariance Matrices
While in this work can be estimated using arbitrary positive definite covariance matrices Σ E Rdxd Such a covariance matrix can be set according the users preference, encoded using the fact that the eigenvectors of Σ represent directions of interest and their corresponding eigenvalues represent levels of smoothing.
isotropic covariance matrices of the form Σ = = -2I,0 0, SmoothHess
In the simplest case the user already has an orthornomal eigenvector basis in mind, (1,...,, Og E Rdd as well as desired levels of smoothing O1....,o Fa > 0. In this case Σ may be set by constructing eigenvector and eigenvalue matrices Q = [.i]...|oo Jaa] E dddd and A = diag(o diagg aa) E gd and simply multiplying = QAQ".
Alternatively, the user may only have k < orthonormal eigenvectors in mind, .1..... Og E Rd for which they wish to smooth at specific levels di...., .k > 0. procedure such as Bram-Schmidi
may be used to find Q4+1,..., 0d E Rd which extends 01 01..... Q2 to an orthonormal basis of Rd Again, the smoothing levels along the eigenvectors d++1 preference. For instance, one may select 0++1 these directions. Just as above may set by constructing eigenvector and eigenvalue matrices Q jail... |aa] E dd and A di..
dg may chosen according to user R 0 if minimal smoothing is desired along
(1.....0.)) da)
It dxd and simply multiplying Σ QAQ".
26

---

E Experiment Setup
E.1
Datasets and Models
In this work we make use of six datasets, two synthetic datasets (Four Quadrants, Nested Interactions) three benchmark datasets (MNIST, FMNIST, CIFAR10) and real world medical dataset (Spirometry). Below we describe these datasets and training details.
Four Quadrant. The Four Quadrant dataset consists of points E IR² sampled from the grid [- x 2- 2] with a spacing of 0.008. A 6-layer fully connected ReLU network was trained using RMSProp [83] on the Four Quadrant dataset achieving a final mean- rrean-squaee-erree of z le-4. Training lasted for 40, 000 iterations with a batch size of 128 and a starting learning rate of 1e-3 which was decayed by a factor of 1e-1 at iterations 5000, 10, 000 and 20, 000.
Nested Interactions. The Nested Interactions dataset consists of points x E R sampled from the grid [-2 2, 2] x 2- with spacing of 01008.. A 6-layer fully connected ReLU network was trained using RMSProp on the Nested Interactions dataset achieving a final mean-squarederroor of R le-1. Training lasted for 200, 000 iterations with batch size of 64 and a starting learning rate of le-3 which was decayed by a factor of le-1 at iterations 40, 000, 80, 000, 120, 000 and 160, 000. For more details on Nested Interactions see App. F.3.
MNIST MNIST consists of 70,000 28x28 greyscale images, each corresponding with one of the digits 0-9. There are 60, 000 and 10, 000 images in the pre-defined train and test sets respectively. We further split the train set into 50, 000 fully connected network with dimensions images for training and 10, 000 for validation. A 5-layer was trained with stochastic gradient descent for 30 epochs, with batch size of 128 and a starting learning rate of 1e-2 which was decayed by a factor of 1e-1 at iterations 4, 000 and 8, 000. Final accuracies of 2 100% (Train) z 100% (Val) and Rs 98% (Test) were achieved. All images are flattened. For PMSE all methods were evaluated on 200 test points. For adversarial attacks all methods were evaluated on 400 test points. The threshold for the quadratic optimization attack was set to T = 0.98.
11--3002200-20--21
FMNIST. FMNIST consists of 70,000 28x28 greyscale images, each corresponding with one of 10 articles of clothing. There are 60, 000 and 10, 000 images the pre-defined train and test sets respectively. We further split the train set into 50, 000 images for training and 10, 000 for validation. Final accuracies of xz 93% (Train) R 93% ( (Val) and 88% (Test) were achieved. The network and training details are identical to that used for MNIST above. All images are flattened. For PM all methods were evaluated on 200 test points. For adversarial attacks all methods were evaluated on 400 test points. The threshold for the quadratic optimization attack was set to T 0.98.
CIFARIO CIFARIO consists of 60,000 3x32x32 RGB color images, each corresponding with an animal or vehicle. There are 50,000 and 10,000 images in the pre-defined train and test sets respectively. We further split the train set into 40, 000 images for training and 10, 000 for validation. A ResNet-18 18 [32] was trained on CIFAR10 for 55 epochs using batch size of 128. The first 5 epochs were used for warmup with a starting learning rate of 1e-2 ending at 0.5. For the rest of training a cosine decay schedule was used [52] decaying down to le-5 by the final epoch. Augmentations used for training were (i) Random Horizontal Flip (p- 0.5) (ii) Color Jitter (p- 0.8) with brightness. contrast and saturation values equal to 0.4 and hue value 0.1 (iii) Random Grayscale (p- 0.2). Final accuracies of & 85% (Train, Augmentations) R 96% (Val, No Augmentations) and R 90% (Test, No Augmentations) were achieved. No augmentations were applied to the validation/test data when evaluating explainers. For ParSE all methods were evaluated on 100 test points. For adversarial attacks all methods were evaluated on 200 test points. The threshold for the quadratic optimization attack was set to T = 0.8
The Spirometry dataset uses raw exhalation curves, measured in volume over time, recorded during a spirometry exam. Each spirometry curve measured in 10ms intervals, which we downsample to 50ms intervals and limit to 15s in total length, resulting in 300 features. We use the UK Biobank dataset, which is a large, population-bas based study conducted in the United Kingdom. Participant statistics have been previously reported in Sudlow et al. [79]. The UK Biobank records 2-3 exhalation efforts for each participant, using a Vitalograph Pneumotrac 6800 device". If two efforts are recorded as passing acceptability criteria and are also reproducible (S 5% difference in Forced Vital Capacity (FVC) and Forced Expiratory Volume in Second (FEV))) then the third
Spirometry.
27

---

Dataser
MNIST
FMNIST
CIFARIO
Function
Class I Logit (
Int. Neuron
Class
Logit
Int. Neuron
Class Logir (1)
Int. Neuron (8)
0.25
0.50
1.00
0.25 0.50 1.00
0.25
0.50
1.00
0.25
0.50 1.00
0.25 0.50
1.00
0.25
0.50
1.00
SH+SG (Us)
4.5e.5
1.8e4
1.8c4 7.2e-4
4.5e:5
1.8e-4
7.2c-4
4.5e.5
188c4 7:2c-4
Llee5
4.6e.5
1.ee4
1.ee5
1.8-4
[76]
4.5e.5
4.8e-4 7:2e
1.8c
7.2e-4
4.5e-5
4.5e
1.8e 7.2e-4
1.le5
4.6e
1.8e
1.le.
4.6e.5
1.8e
G
400.0 200.0 95.0
400.0
200.0
100.0
350.0
160.0
75.0
360.0
190.0
95.0
11.5
6.0
11.0
6.0
4.0
SPG
390.0 200.0
100.0
400.0
200.0
100.0
390.0
180.0
95.0
360.0
190.0
95.0
14.0
8.0
4.0
16.0
7.0
4.0
SW (H+ G
190.0 95.0
55.0
190.0
95.0
50.0
170.0
80.0
45.0
170.0
90.0
45.0
8.5
8.5
8.5
11.0
11.0
11.0
SW G
190.0
95.0
55.0
190.0
95.0
50.0
180.0
100.0
65.0
170.0
90.0
50.0
8.5
8.5
8.5
11.0
11.0
11.0
Table 3: Selected values of o² and 8, achieving the lowest average PMSE on on held-out validation set.
effort is omitted. A common metric used to evaluate lung health is the Forced Expiratory Volume in 1 Second (FEV1), 1), which is the maximum volume of air that can be expelled by the participant in 1 second [59]. Note that participant's FEV, measurement is taken as the maximum FEV1 1 over all recorded exhalation efforts during a single visit.
The spirometer automatically evaluates effort against number of acceptability criteria. One such criteria the detection of coughing. In our experiment, we use a subset of exhalation efforts where coughing was detected and train a CNN predict the participant's final FEV1 measurement. This subset contains 8,721 samples, which we split into training (80%) and test (20%) partitions. We follow the preprocessing in Hill et al. [36] to ensure that participants have at least one effort that passes quality control where FEV, can be measured. The trained CNN includes 10 convolution blocks. Each block contains 1-d convolution of kernel width 200 and 20 channels, batch normalization, dropout (p=0.5), and skip connection. The model is trained using mean squared error (MSE), achieving 0.563 MSE on the train set and 0.547 MSE on the test set.
E.2 Hyperparameters
For ease of reading, 2c, b- c} C R denotes the set of points between a (inclusive) and b (exclusive) at intervals of size c. Here a and b are chosen such that a and (b a) mod An example is: [0.1 1.0: 0.1] = {0.1, 0.2, 0.3, 0.4, 0.6, 0.7, 0.93.
the following notation below: (a b c] fa, + C, 2c,...,b
E.2.1 PMSE
For each of the 18 combinations of dataset, function and neighborhood size the performance of and a are validated on held out set before selection. 200 validation points are used to choose 3 for MNIST and FMNIST and 50 validation points are used for CIFARIO. For SmoothHess + SmoothGrad and SmoothGrad only three values of a are validated, based on the common sense criterion = e/vd: given neighborhood size and dataset with dimension is chosen from a E {=/2vd, 3
The following values of 3 are checked on validation set:
MNIST and FMNIST / E [0.1 1: 0.1] U [1 20 1] U [20 95 : 5) U [100 800 10]
CIFARIO: 8 [0.1 : 1 0.1] [1: : 20 : 0.5] U [20 95 : 5) U [100 800 :10]
The values of a for SH and SH + SG. and 3 for SP H + G and SP G which achieve the lowest PASE on the validation data are shown in Table 3. The results for SoftPlus are interesting: (i) We see that MNIST and FMNIST results. for both class logit and interior neuron, are consistent for fixed e. Further we see that CIFARIO results between class logit and interior neuron are consistent for fixed e. (ii) The optimal value of 3 seems to be approximately proportional to the value of e. While our results in Table show that SmoothHess is better at capturing local interactions than the SoftPlus Hessian, the results in Table 3 indicate that the relationship between f. fo and B warrants further exploration.
E.2.2 Adversarial Attacks
200 validation points were used to choose 3 and の² for MNIST and FMNIST and 50 validation points were used to choose 8 and o² for CIFAR10.
28

---

Dutuset
MNIST
FMNIST
CIFARIO
Atteck M
0.25
0.50
0.75
1.25
1.75
0.25
0.50
0.75
1.25
1.75
0.2
0.3
0.4
SH-SG (Us)
lo.3
Se-3
36-2
Se-2
1.5e-1
lo-3
4e-2
3e-2
8e-2
3e-2
5o-5
lo44
4c-4
4o-4
7.5c-4
SG [76]
3e-4
SPG
11
11
Table 4: Selected values of o² and 3, achieving the lowest post-hoc accuracy of adversarial attacks held-out validation set.
The following values of o² were validated for adversarial attacks:
MNIST and FMNIST: o² {0.001, 0.005) [0.01: 0.1 0.01] [0.15 1.00 0.05)
CIFARIO: d² E (5c-05, 7. 7.5c-05) (/.0000 0.0006
0.006 0.001]
The following values of 3 were validated for adversarial attacks:
MNIST and FMNIST: 3 E [1 :0 U [20 100 5] [100 210 10).
CIFARIO: 3 e [1 10 : 1] U [10: 45
The values of o² and / that achieve lowest validation post-hoc accuracy are shown in Table 4. generally the case that parameters corresponding with increased smoothing (larger o² and smaller 3) achieve better results (lower post-hoc accuracy) for large e, and parameters corresponding to less smoothing (smaller o² and larger (3) achieve better results for small
Following intuition,
Additional Experiments
PASSE
Comparison with Swish: In Table PMSE results are presented for five methods: the first (SmoothGrad) and second (SmoothHess: SmoothGrrd) order Taylor expansions the ReLU network convolved with Gaussian, the first and second order Taylor expansions of the SoftPlus smoothed network and the vanilla (unsmoothed) Gradient. We present Table 5, version of Table which includes additional results comparing with Swish [65] smoothed networks. Swish, an alternative smooth activation to SoftPlus, is formally defined as Sw)(r)
sigm)id(8z:) where sigmoid(z) = T+exp (-x] and 3 is a hyperparamter determining the level of smoothing. It can be seen in Table that Swish is generally less effective then SoftPlus, and is outperformed by our method at each combination of dataset and locality.
We report the standard deviation of the PMSE for each method, dataset, function and neighborhood size E in Table 6. Of the 18 PMSE results, SmoothHess + SmoothGrad attains the lowest standard deviation for 15 and ties with SoftPlus Hessian SoftPlus Gradient for 2. The standard deviation of SoftPlus Hessian SoftPlus Gradient is the lowest for FMNIST internal neuron at = = 0.5, achieving 2.4e-7 while SmoothHess + SmoothGrad achieves 2.5e-7.
Standard Deviation:
We repeat PMSE experiment for the predicted class logits of CIFAR10 using ResNet101, reporting results in Table 7. It can be seen in the leftmost column that our method, SH + SG, achieves superior performance the competing methods at each locality, indicating that SmoothHess can generalize to larger network architectures. The standard deviation of PMSE, as well as choices of o2 / B (as selected from a validation set) are reported in the center and rigthmost columns, respectively.
ResNet101:
F.2 Adversarial Attacks
Comparison with Vanilla Hessian: One may use the vanilla Hessian of the predicted SoftMax probability, which admits higher order derivatives, construct adversarial attacks. Table is an updated version of Table 2 which includes results for attacks using the vanilla Hessian vanilla Gradient (H G), in the third last row.
We see that inclusion the vanilla Hessian generally results in more effective attacks then use of the vanilla gradient alone. H G ties SH+SG and SP H + G for lowest post-hoc accuracy at the
29

---

Dataser
MNIST
FMNIST
CIFARIO
Funetion
Class I Logit (1)
Int. Neuron
Class I Logit
Int. Neuron (L)
Class Logir (1)
Int. Neuron (8)
0.25
0.50 1.00
0.25 0.50 1.00 0.25 0.50 1.00
0.25 0.50 1.00 0.25 0.50 1.00
0.25 0.50 1.00
SH++G (Us)
9.6e-7 7.8-6 6.7e-5
4.9e-8 4.0c-7 3.3e-6 6.5e-7 4.0e-6 4.3e.5
2.0c-8 18c-7 16c-6
9.8e-4 2.2e-2 1.2e-1
8.1e.7 14e-5 16e-4
[76]
4.5e.6
4.1e:5 3.9e
2.1e:7
1.7e.6 1.5e-5 3.0e
2.7e-5
1.0e.7 9.0e:7 7.0e-6
1.3e:2 8.6e: 4.9e-1
1.3e:5 1.le.4 8.3e
SP (H + G)
9.66-6
8.16-5
5.5c-8
4.4c-7 3.70-6
9.6e-7 7.50-6 6.5o-5
3.0c-8
2.1e-7 1.8c-6
2.1c-3 3.3e-2 2.5c-1
1.1e-5
1.0c-4 7.004
SP G
4.6e-64.1e-5 3.9e-4
2.1e-7
1.7c-6 1.50-5 3.2e-6 2.8e-5 2.6c-4
1.0c-7 8.5e-7 7.20-6
L.3e-2 9.0e-2 5.2c-1
5.1e-5
2.9c-4 1.60-3
SW (H+ G)
2.e-6
2.0e-
1.9e-4
1.0e-7
8.3e/7 7.3e-6 1e16 1.7e-5 1.8e-4
5.0e.8 -.0e--843ee7733.. 3,7e6-
I.le22 3.3e-1 7:9e0
6.2e-55
1.7e-3 3.Re-2
SW G
5.6e-6 5.0e-5 4.9e-4
2.4e-7
2.0k-6 1.8e-5 3.96-6
3.5e-5 3.5e-4
1.1e-7 8.6-779.66-8.33----
4.9e-2 9.8e-2 6.0e-1
5.6c-5
3:4c-4 2.0e-
G1731
4.2e-3 1.7e-2 6.7e-2
2.7- 7.0k-3 2.9e-2 3.8e-3 1.5e-2 6.0e-2
1.0c-4 4.0c-4 1.8e-3 3.0e-1 1.2e-0 5.0e-
9.0c 3.5e-3 1.4e-2
Table 5: Average Parsæ results at three radii e. with the inclusion of vanilla Gradient (G). Other methods
include: SmoothHess + SmoothGrad (SH+SG,Us) SmoothGrad (SG) SoftPlus Grad (SP G) SoftPlus Hessian + Gradient (SP + G). Results are provided for the predicted class logit, and the penultimate neuron maximally activated by the "three," dress and cat classes for MNIST, FMNIST and CIFARI0 10 respectively.
Dataset
MNIST
FMNIST
CIFARIO
Function
Class I Logit (1)
Int. Neuron (L)
Class Logit
Int. Nouron (1)
Class Logir (L)
Int. Neuron (L)
0,50 1.00
0.50 1.00 0.25 0.50 1.00
0:25 0.50 1,00 0.25 0.50 1.00
0.25
1.00
SH+SG (Us)
9.9--7 7.3-6 6.4e-5
3.96-8 2.7c77 2.1e-6 L0e-6 7.7e-6 7.7e-6 6.4e-5
2.5e7
2.8e-8
1.9e-6
4.7e-3 I.fe-1 4.3e-1
1.5e-6
3.6e-5 4.2c-4
SG [76]
7.0c-6 7:4-5 7.16-5
2.26-7
1.60-6 1.20-5
5.1e--
4.6o-4
2.0c-7
1.8c-6 1.3o-5
4.76-2 1.9c-1 7.5o-1
4.20-5
4.1c-4 2,96-3
SP (II G)
9.5e
7.7e 5
4.4e 8 2.9e:7 2.1e-6
2.4e
6 5.5e: 1.2e
3.
2.4e-7 1.9e-6 6.5e: 1.3e:1 6.4e:1
5.7e:5
4.0e 3.9e
SP G
7. 7.2e
2.2e:7
1.66:6 1.2e-5
5.4e:5 5.4e:655.e55
6.1e
4.7e
2.0e:7
6.8e-6 4e-5
4.8e:2 1.9e-1 7.6e-1
2.6e
8.9e 3.6e-3
SW (H+ G)
2.6c-6 2.0e-5 1.9e-4
6.5c-8
4.6c-7 4.2o-6
7.2e-6
3.3e-5 3.7c-4
5,3c-8 4.5e-7 4.0e-6
5:0e-2
1.7e0
4.3e1
3.6c-1 1.2e-2 2,7e-1
SW G
1.3e-6 7.2e-5 7.4e4
2.1e-7
1.5e-6
1.2e-5
6.6e-6
6.2e-5 6.2c-4
1.5e-7 1.1e-6 9.0e-6
5.0e=2
2.1e-1 8.6e-1
2.8e-4
1.0e-3 4.2e-3
G
le 8.4e33 3.40-2
7.9e-5 3.0c-4 1.2e-3 2.7e-3 1.1e-2 4.3e-2
9. 4.0e-4 1.4e-3 2.4c-1 9.6e-1 3.70-0
9.t.0099.0......000
3.77e3 1,5e-2
Table 6: Standard deviation of P three radii e. SmoothHess SmoothGrad (SH4SG,Us) SmoothGrad (SG) SoftPlus Grad (SP G) SoftPlus Hessian + Gradient (SP H + G) Vanilla gradient (G). Results are provided for the predicted class logit, and the penultimate neuron maximally activated by the "three," " dress and cat classes for MNIST, FMNIST and CIFARIO respectively.
smallest magnitude (8 = 0.25) for the simplest dataset (MNIST). However, as no smoothing is done, the attacks generated from the H G are generally significantly weaker then those generated using smooth surrogates. In fact, aside from MNIST with E = 0.25, the second order vanilla H + G attacks achieve higher post-hoc accuracy then first-order method SmoothGrad for all datasets and values of e. This is especially apparent for CIFAR1O, the most complex dataset.
F.3 Nested Interactions
We use the Nested Interactions dataset to highlight SmoothHess's ability to capture different interac- tions occurring at various localities around given point. In this experiment we measure interactions around the origin xo- = (0, 0)2 E B...
Just like the Four Quadrant dataset, the Nested Interactions dataset consists of points x E IR2 sampled uniformly from [-2, x 2- 2] c R with spacing of 0.008. We establish different interactions occurring around the origin zo. based upon the distance from To. Specifically, we set the label for a
given point x by: x E (Z.) = y(z)= = t + T2 I E (zo) ( = y(z)=
F1 x E IR I (1200) = y(x) = 222
In words, the interaction between features i1 and 22 is inside the radius-1.2 1.2 ball around zog and is -5 outside of this ball. The interaction between ii and itself is 1 inside the radius-0.6 ball around 250 and O outside of this ball. The interaction between 22 and itself is 0 over all of TR?..
We train 6-layer neural network the Nested Interactions datasct and estimate SmoothHess and SoftPlus Hessian for o² E {1e-6...., lel] and 3 E {le-1...., 4.03 respectively. The interaction results for 21 with itself, z1 with x2 and 22 with itself as function of the level of smoothing (a or (3 are reported in Figure 5.
As the target function y(z) is discontinuous, it is not possible for a network to memorize the Nested Interactions dataset. Thus, there very well may be interactions occurring in the network which are not described as above; the interactions we know occur in the data are not a pure "gold-standard". That being said, Figure 5 shows that SmoothHess captures the interactions as we know occur in the data, and the SoftPlus Hessian does not. This suggests that, to large extent, both (i) the network has
30

---

Value
PMSE
PMSE Std
a², 3 used
0.25
0.50
1.00
0.25
0.50
1.00
0.25
0.50
1.00
SH+SG (Us)
1.3e-4
9.5e-4
6.9e-3
1.7e-4
1.1e-3
1.1e-2
1.1e-5
4.6e-6
1.8e-4
SG
5.3e-4
4.3e-3
3.3e-2
1.0e-3
8.3e-3
6.7e-2
1.1e-5
4.6e-6
1.8e-4
SP (H
3.7c-4
2.6e-3
1.7e-2
6.6e-4
2.7e-3
3.1e-2
41.0
16.0
11.5
SPG
6,2c-4
5.0e-3
3.7e-2
1.le-3
8.8e-3
6.9e-2
55.0
22.5
11.5
Swish (H G)
1.5e-3
1.2e-2
4.5e-2
2.5e-3
8.3e-3
2.8e-2
27.5
17.5
1.5
Swish G
8.8e-4
7.0e-3
5.3e-2
1.4e-3
1.1e-2
8.3e-2
60.0
1.0e6
1.0e6
G
3.8e0
1.5e+1
6.0e+1
3.3e0
1.3e1
5,2e1
n/a
n/a
n/a
Table 7: Using a ResNet101 trained on CIEARIO, the average PasE achieved by SmoothHess + SmoothGrad
(SH SG, Us), SmoothGrad (SG), SoftPlus Hessian + SoftPlus Gradient (SP (H+G)), SoftPlus Gradient (SP
G). Swish Hessian + Swish Gradient (Swish (H+G)), Swish Gradient (Swish G) and Vanilla Gradient (G), evaluated as proxy for explainer quality, Results are reported three radii for the predicted class logit. Left: PMSE results are reported. The lowest value in each column is bolded. Middle: The standard deviation of PMSE is reported. The lowest value in each column is holded. Right: The smoothing hyperparameter (a² for SmoothGrad and SmoothHess SmoothGrad, for SoftPlus and Swish) used is reported. Our method, SH SG. achieves the lowest PMSE: e each the three radii e. This indicates that the interactions SmoothHeess captures improve the model network behaviour, even for large networks such as ResNet101.
Dataset
MNIST
EMNIST
CIFARIO
Mognitude
0.75
0.50
0.75
0.3
1.0
(Uss
46.8
26.5
4.5
[76]
49.3
SP (H1 G)
93.0
81.8
51.5
15.8
7.5
79.8
51.0
27.5
0.8
64.5
42.0
31.0
23.5
7.5
SP
93.3
82.3
53.8
16.3
5.0
79.8
51.5
29.5
7.8
1.0
66.5
47.5
36.0
29.5
8.5
H
93.0
81.8
55.3
19.0
11.8
80.0
50.0
30.3
9.5
2.0
68.0
51.5
40.5
32.5
22.0
[733]
93.3
82.8
56.0
18.5
8.8
80.3
52.3
31.8
11.0
2.5
69.0
51.5
41.0
34.0
21.5
Random
99.8
99.5
99.0
99.0
98.8
99.3
98.0
97.3
95.5
93.8
100.0
99.5
99.0
98.5
96.5
Table Post-hoc accuracy of adversarial attacks performed the predicted SoftMax probability. at five attack magnitudes e, with the inclusion methods include: SmoothHess SmoothGrad (SH + SG, Ours), SmoothGrad (SG), SoftPlus Gradient (SP G), SoftPlus Hessian SoftPlus Gradient (SP (H G)) and vanilla Gradient (G). First order attack vectors are constructed by scaling the normalized gradient by and subtracting from the input. Second order attack vectors found by minimizing the corresponding second-order Taylor expansions.
vanilla Hessian vanilla Gradient + G). Lower is better. Other
memorized the data and (ii) SmoothHess captures the network behaviour while the SoftPlus Hessian does not.
F.4 Qualitative Comparison
We present a visual comparison of the interactions found by SmoothHess with those from other methods. Namely, we consider methods that can be interpreted as the quadratic term in second-order Taylor expansion around smooth surrogate network: SoftPlus Hessian and Swish Hessian.
We show interactions found between super-pin super-pixels of CIFAR10 test images. To this end, we utilize the Simple Linear Iterative Clustering (SLIC) [ algorithm to segment the image into 20-25 super-pixels We sum interactions between each pair of features in each pair of super-pixels, before visualization.
Results are shown in Figure 6 for the predicted class logit of a ResNet18 model for three CIFAR10 test images. Here, each row corresponds to a separate image. Test images are visualized column 1. Columns 2-4 correspond to the three methods. For each image, interactions between one chosen
31

---

Figure 5: Three estimated Hessian elements at zo = (0. 0) for a 6-layer ReLU Network f: R² - trained on the Nested Interactions dataset. (a) SmoothHess (SH) is estimated with isotropic covariance Σ = o²1 using granularly sampled o² {1e-6.... ..... 10). At minute log10 0² -4 either hyper-local noisy behavior captured, or smoothing so negligible that the smoothed function is approximately piece- wise linear with a low-magnitude Hessian, For o² ranging from log =- to log 0² 0 we see SmoothHess reflects the interactions in the dataset: Starting at logı0 o² -4 both .y and 0177 have an interaction ss until the interaction between 21371 hegins to dip to 0 around log,o o² 115.. Finally around log1a a² the interaction 2122 begins to dip toward -5, until log o² 0 when o² is large that samples outside SmoothHess estimation. (b) The Hessian of the SoftPlus SH V2U** (2)(0) vs log
SP Hess
V2ff(0) vs logio B
(a)
(b)
Smoother
SoftPlus Hessian Values
Smoother
SmoothHess values
8
6
3.
450
120
02
10g
10g /
training data distribution are incorporated
smoothed function fo (SP Hess) is computed using granularly sampled / E (le-l..... 1e4]. Here, as is decreased, it is not apparent that the variety of interactions in the Nested Interactions datasct are captured, cither in relative ordering magnitude.
super-pixel (outlined in black) and each other super- super-pixel are visualized as heatmap overlaid upon the image. In order to facilitate comparison across images and methods, standardize the heatmap colorbar to range between the most negative and most positive interaction values on per-image and method basis.
One interesting trend seen in each case is that there a strong positive interaction between the chosen super-pixel and one other super-pixel which (a) is spatially nearby and (b) contains the class object of interest. For example, front wheel. In the second row, the tip of the frogs head can be seen interact positively with the side of the head. In the third row, the upper and lower portions of the dogs front leg have a strong positive interaction.
the first row, the side-view mirror of the car positively interacts with the
Due to the subjectivity this comparison, we include quantitative results above each image. Specif- ically, we indicate the PMSE each method achieves within an 0.25 ball around each image. Optimal smoothing parameters were chosen for each method for this task (see Table 3). can seen that SmoothHess achieves the lowest PASE in each case by wide margin. Thus, SmoothHess may be the preferable option one wishes for visualization which best reflects the network's behaviour in an 0.25 ball around the image.
32

---

Test Image
SmoothHess
SoftPlus Hessian
Swish Hessian
S MAE 2.0e-4
8 UsE 7.0c-4
S MSE 1.5e-3
Max Poottive Interaction
MSE 1.00-44
P usE 3.0e-4
9 MsR 7.0e-4
No Interaction
5 urE = 1.0e-4
S usE = 4.0e-4
MSE 8.0e-4
Nogative
Figure 6: Visualization of interactions between super-pixels found for ResNet18 on CIFARIO by SmoothHess. SoftPlus Hessian and Swish Hessian. Results shown for test images of frog and dog in the first second and third rows respectively. Each image is visualised in column one. Images are segmented into 20-25 super super-pixels using the SLIC algorithm [2]. Interactions are summed between each pair of features in each pair of super-pixels. We show interactions with one given super-pixel in cach image, outlined in black. SmoothHess, SoftPlus Hcssian and Swish Hessian interactions for this super-pixel are visualized heatmaps overlaid upon the image in columns two, three and four, respectively. The heatmap colorbar standardized range between the minimum comparison. Quantitative PASE results for e 0.25 are shown above each method-image pairing, with Hesothlesss achieving the lowest PASE in all cases.
maximum interactions on each image- method pair separately, facilitate
33


Reviews:
[{"review_id": "doqjs1dy02", "reviewer": "Reviewer_3HE6", "summary": "The author introduces a novel method to approximate the Hessian of a smoothed\nloss function at a specific point. By incorporating a variation of Stein’s lemma,\nthey cleverly estimate the gradient over a smoothed surrogate of a ReLU network\nat that point. To gauge the method’s effectiveness, they measure its performance\nby testing the normalized L2 distance between the loss function and its first-\norder and second-order approximations within a ball centered around the given\npoint.", "questions": "1. line 103 \"An general L-hidden...\" is a typo?\n\n2. line 287 \"with the non-smooothed...\" is a typo?", "limitations": "The article address the time and space complexity of the method in the presence of a high-dimensional input.", "rating": 7, "confidence": 4, "soundness": 4, "presentation": 3, "contribution": 3, "strengths": "The method offers a straightforward and comprehensible approach, presenting\na powerful and novel tool for approximating the Hessian of a ReLU network at\na specific point. The author’s analysis is comprehensive and remarkably easy\nto follow. The experiments conducted in the paper appear relevant, and the\nobtained results show promising potential.", "weaknesses": "While the ideas and methods of the paper are solid, I came across some read-\nability issues.\n\n1.\nThe statement of Theorem 3 (line 211) is confusing. The real positive\n$\\delta$ is used without declaration, while the random variables $\\delta_1,\\ldots,\\delta_n$ are declared. This makes the reader wonder how one can compare a probability with a random vector, and what is the meaning of this inequality which is essentially a random variable itself. I had to read the proof in the appendix to fully understand the statement of the theorem. I think using a different symbol for $\\delta$ might increase readability.\n\n2.\nI would add the definitions of the oracles used in the article for clarity.", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "\nWe appreciate your positive assessment of our work. Below, we address the weaknesses and questions that you have brought up. \n\n\n## Weaknesses\n\n---\n\n**W1: Regarding the statement of Theorem 3.**\n\nThank you for your careful reading of the statement of Theorem 3. We agree with the reviewer’s suggestion. In the camera ready version, we will change the real positive $\\delta$ to $\\gamma$ in the statement and proof of Theorem 3. We will also explicitly declare $\\gamma$ in the statement of Theorem 3. \n\n**W2: Regarding oracle definitions**\n\nA: We will add the definitions of both oracles to the camera-ready version. Namely, a zero-th order oracle is a function call. That is, it is an oracle that, given a function $f : \\mathbb{R}^d \\rightarrow \\mathbb{R}$ and input $x \\in \\mathbb{R}^d$, returns $f(x)$. In the context of our work, this amounts to a neural network forward pass. Likewise, a first-order oracle is an oracle which, given $x$, returns $\\nabla_x f(x)$. In the context of our work, this amounts to a network backpropagation call. \n\n\n## Questions\n\n---\n\n**Q1 and Q2: Typos: “An general L-hidden … “, “non-smoothed”**\n\nThank you for catching these. We will change the sentence to “A general L-hidden … “ and “non-smoothed” to “unsmoothed” in the camera ready version.\n"}, {"review_id": "pUGmqvR13j", "reviewer": "Reviewer_FahQ", "summary": "The work proposes an approach for measuring feature interactions for neural networks that have no\nuseful second order gradients or Hessians (notably those employing ReLU). The approach convolves a\nsmooth normal with the network as a proxy with usable Hessians and then offers a means of\nestimating them using sampling of first-order gradients, thereby making it more tractable than\nmethods that rely on computing Hessians naively. Experiments evaluate the proposed method and\nseveral adaptations of existing methods for computing 1) the 3rd-term (second order gradient) of\nthe Taylor expansion of the output of a neural network around an input, and 2) for adversarial\ninputs via optimization using (1). In these experimental terms, the proposed approach (aided by\nSmoothGrad for the 2nd Taylor term) performs better than those based on only the first order\ngradient or those that use include both first and second order gradient but proxy ReLU models using\nSmoothPlus activation.", "questions": "- Conceptually, when sampling as per MC-estimation, if samples end up not spanning more than 1\n  linear region of the original network, would the estimate be zero in the same sense as second\n  order gradient is zero?\n\n- What is the impact of smoothness/covariance on sample bound as per Equation 9?\n\n- Does \"SP H + G\" mean \"SoftPlus (Hessian + Gradient)\" or \"(SoftPlus Hessian) + (ReLU Gradient)\"?\n  Line 285 suggests the latter but the label/caption suggest the first. If the first, the\n  comparison may not be as fair as it could be.\n\n- Example of Figure 3 is difficult to understand. Why are the arrows labeled \"interaction\" crossing\n  from one method to another? Also, SmoothGrad vs. SmoothHess is indicated by both color and the\n  separate graphs? Did you mean that green/red are indicators of positive or negative interactions\n  for both of those methods or for one respectively? Or perhaps the second of each pair of graphs\n  shows the interaction of the other timesteps with the first timestemp only? Also, what is the\n  kernel width for the example shown in Figure 3?\n\n- The covariance matrix is set in a specific way in the text but the benefit of the option to set\n  it otherwise are mentioned. How would one go about setting it alternatively and for what\n  purposes?", "limitations": "The work describes experiments with adversarial attacks which could use some mention with relation\nto negative societal impact.", "rating": 6, "confidence": 2, "soundness": 3, "presentation": 2, "contribution": 3, "strengths": "+ Almost strictly better experimental results than the compared-to methods.\n\n+ Hessians for networks which don't naively have them due to technicalities should be generally\n  useful in a wide variety of applications.", "weaknesses": "- Experimentation focuses on agreement between a network and one evaluated as per Taylor expansion\n  using SmoothHess. While this may be a fairly general test, it does not offer any insight into\n  benefits in cases where the method would actually be useful to use, with no trivial alternatives.\n  Presumably in most cases we can already evaluate the network regardless of Hessian non-existence\n  so at surface level there is no need for Taylor expansion and thus the Hessian of the smooth\n  proxy. The adversarial attack is also based on use of Taylor expansion of smoothed model as proxy\n  so it does not improve this problem. The final case-study does present interaction results but is\n  difficult to follow (see questions below). I suggest, then, that the paper include evaluation of\n  the methods that focus on the Hessian, instead of network output. One example is to test\n  adversarial attacks on gradient-based attributions by gradient descent on the gradient (i.e.\n  using Hessian). Compare to other methods for attacking attributions.\n\n- Tradeoffs with respect to the parameters are not explored or described. The use of smoothing\n  serves to increase the impact of network regions further away from each point but smoothing is a\n  double-edged sword. If smoothing (and subsequent sampling) is insufficient compared to the sizes\n  of linear network regions, the results may not be useful. On the other hand, large enough\n  smoothing reduces the network to a constant.\n\n  One impact of this tradeoff is how fast gradient descent can arrive at good solutions. Thus in\n  addition to the above suggestion of testing usefulness of the Hessian, I also suggest that the\n  effect of the parameters (or smoothness) have on that usefulness. For example, sufficiently\n  significant smoothness I expect to result in useless gradients of gradients for attacking\n  attributions.\n\nSmaller things:\n\n- What is δ (without subscript) in Equation 8a? Did you mean δᵢ ?\n\n- On line 117, the expectation and perturbation are presumed not significant to the point that the\n  Hessian is 0 almost everywhere.\n\n- By \"second-order interactions\" did you mean interactions via second-order gradients? Interactions\n  seem already second-order so saying second-order interactions is suggesting third-order\n  gradients.\n\n- In Table 1, are the two best (bolded) results in the 4th numerical column identical or is it just\n  that the printout identical?\n\n- Around line 255, you write \"We denote Δ* ... as the output after x₀ is attacked\". This may be a\n  bit confusing as I think you mean Δ* as the perturbed input, not the network output. Is that\n  right?\n\n- Line 349 writes \"Interestingly, the smaller kernel width constraints the interactions to features\n  within a small locality.\" Why is that interesting? Isn't that expected?\n\n- Typo near \"between between\".\n\n- The Adversarial Attacks paragraph on line 254 and Equation 12a is missing the target class of the\n  attack.\n", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "Thank you for reviewing our work. **We've addressed your main concerns in our official rebuttal. Other concerns have been addressed in a comment due to the character limit.** \n\n---\n\n**W1  Experimentation Focus** \n\nWe emphasize that, while we agree as to the potential for use in optimizations, SmoothHess is foremost meant to quantify feature interactions. Assessing how well any given interpretability method explains network behavior is non-trivial: it is often subjective and dependent on expert knowledge. A number of techniques have been introduced to evaluate univariate methods[1,2]. Note that these methods don’t necessarily correspond to real world uses. They are used as a quantitative proxy for interpretation quality, but are not the final goal. \n\nAs SmoothHess and SmoothGrad together define a second-order Taylor expansion we developed P-MSE and the adversarial attack experiment as intuitive quantitative measures of interaction quality, using this Taylor expansion as a proxy. Low values of P-MSE may indicate that the average interactions are captured. Successful adversarial attacks may indicate that extremal interactions are captured.\n\nWe find your experiment valuable, but believe it may be subject to the same comments given regarding our adversarial attack experiment. While the use of SmoothHess to attack the function f presupposes access to f, the use of SmoothHess to attack $\\nabla f$ presupposes access to $\\nabla f$, which is used for SmoothHess computation. We believe both experiments are interesting, but note that the adversarial attack on f may be more in line with our goal of interpreting interactions for f. This is as opposed to interpreting first order importance for $\\nabla f$, which is also interesting. We leave this for future work.\n\nAs stated above, the main application of SmoothHess is to be used to interpret network behavior. This is exemplified by the spirometry case study, which we clarify in the global response(see G2).\n \n**W2  Smoothing Tradeoffs**\n\nWe agree with this need and note that we consider the benefits of SmoothHess to be highlighted from this perspective.\n\nFirst, from our viewpoint, the trade-off from smoothing is not inherent, i.e. “smoothing level $\\sigma$ is too high for a given point/network”, but rather predicated on the information one is searching for. E.g., assume one is given a point $x \\in \\mathbb{R}^d$ and network $f : \\mathbb{R}^d \\rightarrow \\mathbb{R}$, and they wish to find the interactions in the radius 0.001 ball: $B_{0.001}(x) \\subset \\mathbb{R}^d$. If there exists a linear region $Q \\subseteq \\mathbb{R}^d$ s.t. $B_{0.001}(x) \\subseteq Q$ then, from the Hessian perspective of interactions, there is _no interaction occurring_ in $B_{0.001}(x)$. In this scenario using a tiny $\\sigma$ would actually be appropriate for answering the users’ question, _precisely because it would result in a near zero Hessian._ \n\nHerein lies the relative benefit of our method over SoftPlus. Given a ball of radius r over which a user wishes to model interactions, $\\sigma$ can be chosen so that SmoothHess draws samples from near the boundary of this ball. This is because the d-dimensional Gaussian with covariance $\\sigma^2 I$ converges to a Uniform distribution over the sphere of radius $\\sigma \\sqrt{d}$ as d goes to $\\infty$ [3]. Thus the covariance may be set to $\\frac{r}{\\sqrt{d}}$, ensuring that SmoothHess uses samples near this sphere. This is explained in lines 171-177. \n\nAlthough this convergence occurs as d goes to $\\infty$ we have observed that it is close to reality with the datasets used, d = 300,784,3072. We validate this approach in our P-MSE experiment: SmoothHess with $\\sigma$ chosen based upon $\\frac{r}{\\sqrt{d}}$ achieves the best P-MSE over competing methods. Alternatively, the smoothing parameter $\\beta$ for SoftPlus has no clear-cut connection to a locality, for exploitation. \n\nLastly, we point to the Nested Interactions Experiment in App. F.3, in which we assess the ability of SmoothHess to pick up interactions as $\\sigma$ is varied. We sample a set of points $x_1, \\ldots, x_N \\in \\mathbb{R}^2 $ uniformly from $[-2,2] \\times [-2,2]$ and set targets $y(x)$ by: $x \\in B_{0.6}(0) \\implies y(x) = \\frac{1}{2} x_1^2 + x_1 x_2$, $x \\in B_{1.2}(0) \\backslash B_{0.6}(0) \\implies y(x) = x_1x_2$  $x \\in \\mathbb{R}^2 \\backslash B_{1.2}(0) \\implies y(x) = -5x_1x_2$. We train a network to “memorize” this dataset.\n\nThus, we have created “ground truth” interactions for _network behavior_ which change with distance to the origin. We estimate SmoothHess at $x_0 = 0$ with $\\sigma \\in \\{1e-6, \\ldots, 10 \\}$ and plot the interactions in Figure 5a. The interaction between $x_1$ and $x_2$ should be $\\approx 0$ within a single region. At slightly larger localities one should expect to capture the noisy behavior of the network. Then an interaction of $1$ should occur and then a decrease to $-5$. Our results in Figure 5 capture this behavior. The behavior at largest $\\sigma$ is due to the many samples incorporated outside of the train distribution. \n\n---\n\n**S4 Precision**\n\nThis is due to precision. Checking with higher precision, we found SmoothHess P-MSE = 4.9e-8 and SoftPlus Hess P-MSE = 5.5e-8. Thus, SmoothHess achieves lowest P-MSE in all cases. We will add this to the caption.\n\n**Q3 Fair comparison**\n\nWe use SP H + SP G. We feel this choice is more fair. Table 1 in the response pdf shows SP G significantly outperforms (ReLU) G. We will make this clear in the final version.\n\n---\n\n**References**\n\n[1] Hooker, Sara, et al. \"A benchmark for interpretability methods in deep neural networks.\" Advances in neural information processing systems 32 (2019).\n\n[2] Chen, Jianbo, et al. \"Learning to explain: An information-theoretic perspective on model interpretation.\" International conference on machine learning. PMLR, 2018\n\n[3] Roman Vershynin. High-dimensional probability: An introduction with applications in data science, volume 47. Cambridge university press, 2018. "}, {"review_id": "vkZOmMOZOH", "reviewer": "Reviewer_GnHq", "summary": "The paper propose to compute Hessian of a ReLU neural network by a convolution of the network function with a Gaussian distribution function. It's equivalent (or cover) some prior work, e.g., SmoothGrad. The method is generally easy to understand, without modifying the network structure. A non-asymptotic bound on the sample complexity is provided.", "questions": "Some notations are confusing, e.g., caption of figure 2, is it $\\log(10\\sigma)$ or $\\log_{10}{\\sigma}$? \n\nalso the text, \"Aside from at minute\", what does \"minute\" mean here?", "limitations": "See above", "rating": 6, "confidence": 2, "soundness": 3, "presentation": 2, "contribution": 2, "strengths": "Strengths:\n-- The proposed smooth-surrogate is easy and straightforward to understand, the computation method is somewhat novel using Stein's Lemma, following [40] with a mild extension to all Lipschitz continuous functions, which addresses the challenge posed by the zero Hessian of piecewise-linear ReLU networks. The paper provides a non-asymptotic bound on the sample complexity of the estimation procedure. \n\n-- Proposed SmoothHess can be applied post-hoc and does not require any modifications to the ReLU network architecture. This allows analysis and interpretation of the interactions of existing ReLU neural networks without retraining or redesigning the model.\n\n-- Efficient sampling: SmoothHess utilizes an efficient sampling algorithm, which only requires network gradient calls, making it computationally feasible to large-scale networks and complex datasets.", "weaknesses": "Weaknesses:\n-- Scope limited to ReLU networks: While SmoothHess addresses the challenge of zero Hessian in ReLU networks, its applicability is limited to this specific activation function. It would be useful to discuss the potential generalizability of the method to other types of neural networks with different activation functions.\n\n-- Generalize beyond simple datasets: SmoothHess can be applied to large models and datasets with its simple computations, it would be beneficial to clarify the generalizability of SmoothHess and the potential limitations in terms of dataset diversity.", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "Thank you for taking the time to review our work. Below we respond to the weaknesses and questions brought up. \n\n---\n\n## Weaknesses\n\n**W1: Limitation of scope to ReLU**\n\nWe appreciate the comment as we feel this is an important point to clarify in our work. Our method is not limited only to ReLU networks and in fact it is of use in far broader settings. We focus our presentation on ReLU networks because, for their logits and internal neurons, smoothing _is required_ for the computation of non-zero Hessians. Proposition 1 states an equivalency between the expectation our method estimates and SmoothHess for _any Lipschitz continuous function_. Thus our method is applicable to any Lipschitz continuous network or model for which one has access to the gradients (i.e. can backpropagate through). This covers a wide range of networks, for instance networks which may be composed of convolutional layers, fully connected layers, max pooling, common activation functions such as sigmoid, skip connections and batch normalization among other frequently used network building blocks [1].\n\nSmoothHess still provides useful information for Lipschitz networks which have a non-zero Vanilla (unsmoothed) Hessian. Although the Vanilla Hessian does provide information pertaining to the behavior of such functions, the locality of this information is infinitesimally small. On the other hand, SmoothHess has a smoothing parameter $\\Sigma$ which may be chosen to model the interactions over neighborhoods of varying size, providing the user great flexibility. \n\nWe plan to emphasize this point in our camera-ready paper. \n\n**W2:  Generalization to larger datasets and models**\n\nIndeed, the input size matters, which we point out in the limitations section. Specifically, we note that the most computationally expensive aspect of SmoothHess is the outer-products, which are $\\mathcal{O}(d^2)$. For large $d$, this becomes expensive. We gave the example of ImageNet, which has $d^2 \\approx 10^{10}$. However, we do not view this as a limitation relative to other works, since interaction methods are generally $\\Omega(d^2)$ [2,3,4], as each pair of $d$ features must be compared. In our work, we have run SmoothHess for inputs of up to size $d=3072$, although we have not exhaustively checked the limit for which SmoothHess is feasible. From a qualitative perspective of dataset diversity, we have managed to capture feature interactions for a number of different datasets: MNIST, FMNIST, CIFAR10, real world Spirometry regression, synthetic Four Quadrant and the synthetic Nested Interactions dataset found in Appendix F.3. \n\nWe also wish to highlight the new experimental results reported in the global response (G4), which provide evidence that SmoothHess can generalize to larger models. Specifically, our results show that SmoothHess can capture feature interactions for a ResNet101. Table 2 of the global response pdf shows our results for this experiment. \n\n---\n\n## Questions \n\n**Q1: Figure 2 caption.**\n\nThank you for pointing this out. We mean $log_{10} \\sigma$ and will add the subscript to make this more explicit in the camera ready version. \n\n**Q2: Meaning of “minute”**\n\nHere minute means “very small”.\n\n---\n\n## References \n\n[1] Gouk, Henry, et al. \"Regularisation of neural networks by enforcing lipschitz continuity.\" Machine Learning 110 (2021): 393-416.\n\n[2] Joseph D Janizek, Pascal Sturmfels, and Su-In Lee. Explaining explanations: Axiomatic feature interactions for deep networks. fJ. Mach. Learn. Res., 22:104–1, 2021. \n\n[3] Michael Tsang, Sirisha Rambhatla, and Yan Liu. How does this interaction affect me? inter- pretable attribution for feature interactions. Advances in neural information processing systems, 33:6147–6159, 2020. \n\n[4] Dhamdhere, Kedar, Ashish Agarwal, and Mukund Sundararajan. \"The shapley taylor interaction index.\" arXiv preprint arXiv:1902.05622 (2019).\n"}, {"review_id": "iel30N6s6B", "reviewer": "Reviewer_HbQW", "summary": "This paper introduces SmoothHess, a method to define a smooth approximate Hessian for networks where many components of the Hessian may be zero (e.g. ReLU networks). This can aid in tasks where interactions between features are needed, but where those components of the exact Hessian are zero. SmoothHess is defined as the convolution of the output of the network with the Gaussian. It is easy to show, using Stein's lemma, that this is equivalent to the expectation of the gradient of a perturbed output over that Gaussian distribution. They then do a number of experiments. In the first two sets of experiments (perturbation MSE and adversatrial attakcs) on three datasets to show the advantages of using SmoothHess over alternative smoothings of ReLU, such as softPlus. In both experiments, SmoothHess significantly outperforms SofPlus. They also test this on a 4 quadrants dataset with similar results. They also apply this to a real-world medical dataset of spirometry to measure feature interactions.\n\n", "questions": "1. Why only consider diagonal Hessian components, while off-diagonal blocks of ReLU networks can be non-zero?\n2. How do baselines other than SoftPlus (e.g. swish) perform?\n3. In the spirometry experiment, what is the conclusion? Can you trust or validate that the interactions revealed by SmoothHess are real? Does it help predict which datapoints to discard, or is it about which ones not to drop? \n4. regarding the quadratic computation of SmoothHess, is there an alternative for feature interactions which is not quadratic? \n", "limitations": "The paper does discuss limitations, notably the quadratic cost in feature dimensions, which is prohibitively expensive for large models/datasets. I don't see an immediate societal impact issue.  ", "rating": 6, "confidence": 3, "soundness": 3, "presentation": 4, "contribution": 3, "strengths": "The paper shows the benefit of using a smoothened Hessian of ReLU networks for tasks involving feature interactions. The benefit of the method is that it is post-hoc and can be applied to the output of a trained network without a modification to the architecture. The theoretical derivations are easy to follow and well organized. The experiments also make a good case for using their construction, showing superiority to SoftPlus.\nOverall, the simplicity of the method is its main strength.", "weaknesses": "My main reservation is how impactful the work can be. It is a neat idea, but with a rather limited scope. For example, given the expense of computing SmoothHess and the fact that we only need it for diagonal blocks of the Hessian of ReLU networks (explained next) makes the impact rather limited for me. Maybe the authors can explain or emphasize more the breadth of the problems where this would be useful or required?   \nOne minor issue to raise is that, as I understand it, the Hessian of ReLU networks is *not* zero almost everywhere. In multi-layer networks, only the Hessian components involving weights in the same layer are zero. All cross layer Hessian terms should be nonzero. Let $W^l_i$ denote a component indexed by $i$ in the weights of layer $l$. The Hessian term   \n$${\\partial^2 F(x) \\over \\partial W^l_i \\partial W^k_j} $$   \nis in general not zero, even in ReLU networks. \nIn other words, the paper's claims about vanishing Hessian only hold for block diagonal components of the Hessian. They need to argue for why these components are more important than off-diagonal blocks.  \nAdditionally, the baselines in the experiment are rather limited. SoftPlus is not the only smooth alternative to ReLU. The swish activation and similar functions could have been compared as well. \n", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "Thank you for taking the time to review our work. Below we respond to the weaknesses and questions brought up. \n\n---\n\n## Weaknesses\n\n**W1: Impactfulness and non-zero Hessian**\n\n\nWe kindly wish to clarify any misunderstandings about the Hessian we are computing and its use in our work. We are interested in the Hessian used in interpretability [1,2,3,4] which is the Hessian with respect to the input $x \\in \\mathbb{R}^d$, with elements of the form $\\frac{\\partial^2 F(x)}{\\partial x_i \\partial x_j}, \\ i,j \\in \\{1, \\ldots, d\\}$. The interpretation of this quantity is that, following the intuitions from calculus, $\\frac{\\partial^2 F(x)}{\\partial x_i \\partial x_j}$ represents the interaction between features $x_i$ and $x_j$ in affecting the models prediction for any specific point $x$, i.e. $F(x)$. \n\n\nIn the particular comment the reviewer refers to the Hessian with respect to the weights ($\\frac{\\partial^2 F(x)}{\\partial W_i^l \\partial W_j^k}$), which is also a very interesting quantity  and important for optimization, e.g. second-order methods. However in this work we are not interested in this particular Hessian. In this respect, we neither need to assume nor have a block structure. \n\n\nFor the fact that the Hessian (with respect to inputs) of ReLU networks is 0 almost everywhere, please see, e.g. [1]. We convey a proof sketch here: First observe that any affine function $f : \\mathbb{R}^d \\rightarrow \\mathbb{R}^c$ of the form $f(x)  = Wx + b, \\ x \\in \\mathbb{R}^d, W \\in \\mathbb{R}^{c \\times d}, b \\in \\mathbb{R}^{c}$ has a Hessian that is 0. From here one may see that any piecewise linear function has a 0 Hessian, wherever it is defined. Further, note that the composition of piecewise linear functions is piecewise linear. Deep ReLU networks are proven to be compositions of piecewise linear functions [5].  \n\n**W2 : Limited scope of baselines**\n\nFollowing the reviewer’s suggestion we ran new perturbation mean-squared-error (P-MSE) experiments comparing with Swish smoothed networks. We kindly point you toward the global response (G3) for an explanation of our experiment and Table 1 of the global response pdf for the results. SmoothHess outperforms Swish in all 9 permutations of the P-MSE experiment we have run. \n\n---\n\n## Questions \n\n**Q1 and Q2**\n\nPlease see the responses to W1 and W2 above, respectively. \n\n**Q3: Regarding Spirometry Case-Study** \n\nWe kindly point you toward the global response (G2) for clarification on the spirometry experiment, and the conclusions that may be drawn from it. \n\n**Q4: Regarding quadratic cost** \n\nMethods which aim to quantify second-order interactions between every pair of features for $d$-dimensional inputs are generally $\\Omega(d^2)$ given each pair of $d$ elements are compared. We are not aware of works which circumvent this inherent issue. \n\n---\n\n## References\n\n[1] Joseph D Janizek, Pascal Sturmfels, and Su-In Lee. Explaining explanations: Axiomatic feature interactions for deep networks. J. Mach. Learn. Res., 22:104–1, 2021. \n\n[2] Michael Tsang, Dehua Cheng, Hanpeng Liu, Xue Feng, Eric Zhou, and Yan Liu. Feature interaction interpretability: A case for explaining ad-recommendation systems via neural interaction detection. In International Conference on Learning Representations, 2019. \n\n[3]Samuel Lerman, Charles Venuto, Henry Kautz, and Chenliang Xu. Explaining local, global, and higher-order interactions in deep learning. In Proceedings of the IEEE/CVF International Conference on Computer Vision, pages 1224–1233, 2021. \n\n[4] Sahil Singla, Eric Wallace, Shi Feng, and Soheil Feizi. Understanding impacts of high-order loss approximations and features in deep learning interpretation. In International Conference on Machine Learning, pages 5848–5856. PMLR, 2019. \n\n[5] Randall Balestriero and Richard Baraniuk. Mad max: Affine spline insights into deep learning, 2018.  \n"}, {"review_id": "9gvpuhiZKs", "reviewer": "Reviewer_yoRr", "summary": "The goal of this paper is to model the point-wise interactions between input features in a ReLU network. ReLU networks are piece-wise linear, which inhibits using Hessian to model feature interaction (since it’s zero almost everywhere). This paper proposes SmoothHess, a new method to compute feature interaction by the Hessian of a smooth surrogate of the ReLU network. The authors prove that SmoothHess can be estimated using only gradient calls, and prove non-asymptotic sample complexity bounds.\n\nOne result that’s particularly elegant is the connection between SmoothHess and SmoothGrad: combined, they define a local, second-order Taylor expansion of the network. The author designed several evaluations of the proposed method based on this connection.\n\nThe main baseline that SmoothHess is compared against is SoftPlus Hessian which is based on an alternative smooth surrogate of ReLU. These methods are evaluated as local second-order estimates of the ReLU model: they are evaluated by the accuracy of their local approximation, as well as the capability to guide adversarial attacks. SmoothHess performs best in both tasks on three standard image classification datasets. The authors also performed qualitative analysis on a regression problem from the medical domain.", "questions": "Regarding Equation 10: if my understanding is correct, what Equation 10 does is it defines how SmoothHess (and competing methods) can be used in conjunction with gradient methods to define a local, second-order model of the underlying ReLU model. This equation is then used as the basis of experiments which evaluates the goodness of fit from several aspects. Is this understanding correct? Are there alternatives to Equation 10? Is there an alternative way to do similar experiments without Equation 10? Is it possible for this specific definition (for example by using f instead of g as the first term) to bias the results?", "limitations": "The authors have adequately address the limitations.", "rating": 8, "confidence": 4, "soundness": 3, "presentation": 4, "contribution": 4, "strengths": "Although feature interaction is a well-studied topic, I find the proposed method to be quite novel. The motivation of the problem, i.e., the piece-wise linearity of ReLU networks is very clear and the survey of related work seems thorough.\nThe proof sketches provided in the main paper are explained clearly, with detailed references to pieces of proofs from other work.\nOverall, I find the paper to be very well-written and easy to follow. It strikes a great balance between having detailed theoretical results and clearly conveying the intuition.", "weaknesses": "The main weakness of the paper is the comprehensiveness of experiments. Due to computational cost of computing SmoothHess, all experiments are conducted on relatively small-scale datasets and small-ish models. Although I do not see reasons why this method would not work on larger models, it would make the paper much stronger to have empirical evidence.\n\nI also find the qualitative analysis a bit dry. Part of this is due to my lack of familiarity with the Spirometry regression task. But I cannot confidently evaluate the significance of that result, in particular regarding what it means for interpretability of ReLU networks in general. Can this method improve tasks that can benefit from interpretability, like error detection, trust calibration, debugging, auditing, or knowledge discovery? I'm not sure.\n\nThis is a bit disappointing since interpretability is a primary motivation for this work. It would help to at least have some qualitative comparison (and with more competing methods) on more familiar datasets to at least help us understand the behavioral differences between SmoothHess and existing methods.", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "We appreciate your enthusiasm about, and careful reading of, our work. We are especially glad to note your appreciation of the connection between SmoothHess and SmoothGrad. Below, we respond to the weaknesses and questions you have brought up. \n\n---\n\n## Weaknesses\n\n**W1: Generalization to larger models** \n\n\nFollowing your recommendation, we ran an experiment validating that SmoothHess can capture interactions for a ResNet101 trained on CIFAR10. We kindly point you toward the global response (G4) for an explanation and Table 2 of the global response pdf for results.  Thank you for bringing this concern up, this new evidence improves the generalizability of our work.\n\n**W2: Spirometry regression and qualitative comparison**\n\nWe provide an updated explanation of the spirometry case-study in the global response (G2).\n\nWe agree that there are many possible real-world applications for interpretability methods. As our quantitative results in Section 5.2 show, SmoothHess outperforms existing methods of estimating Hessians as well as first-order approaches, suggesting that using SmoothHess should improve downstream applications that use Hessian-based interactions. In addition, SmoothHess is very flexible in that $\\sigma$ may be set to capture local interactions very effectively, at many different localities around a point. Thus, there is some reason to believe SmoothHess may outperform other methods in applications where the Hessian is important; e.g. validating that the model is picking up on interactions which domain experts know to be in the data.\n\nWe agree that having additional qualitative comparisons will be beneficial. We plan to add a qualitative comparison, on image examples, with other methods to the camera-ready Appendix.  Specifically, we plan to segment images, agglomerate SmoothHess interactions within each super-pixel, and carefully compare the difference in behavior with other methods.  \n\n---\n\n## Questions\n\n**Q1: Equation 10. second-order Taylor expansion**\n\nWe appreciate this thoughtful question regarding our experimental methodology. \n\nYour understanding of Equation 10 and its use in our experiments is correct. As you hint at, the most straightforward modification to Equation 10 would be using $g(x_0)$ instead of $f(x_0)$ as the first term. Our motivation in using $f(x_0)$ for the first term stems from the fact that our goal is to evaluate the ability of the smooth surrogates _corresponding gradient-hessian explainer pair_ to model the network locally. We do not consider the quantity $g(x_0)$ (the smoothed function output) to be an explainer of interest; it does not indicate feature interactions such as SmoothHess or feature importance such as SmoothGrad. For this reason we elected to use $f(x_0)$ as the first term for all models. From our perspective this, in some sense, actually removes the bias that could be introduced by the term $g(x_0)$ which we do not aim to evaluate. This allows us to isolate the impact of the gradient-hessian pair in modeling local changes in the function. \n\nThis is opposed to the different, but related, goal of _assessing the smooth surrogates second-order Taylor expansion_  in its ability to model the network, which would require replacing $f(x_0)$ with $g(x_0)$.  \n\nWhile the experiments based on Equation 10 seem the most natural to us, other related experiments can also be run without Equation 10. For instance, as Reviewer FahQ mentioned, one may leverage the interpretation of the Hessian as a first order approximation of the gradient (or some scalarization of the gradient) to create experiments which model the gradient on the first-order, as opposed to modeling the network itself on the second-order as we have focused on. We find Equation 10 to be more in-line with our interpretability goals as generally we hope to quantify the effects of interactions for the function $f$. "}, {"review_id": "1NJUGVgmba", "reviewer": "Reviewer_bZYQ", "summary": "For smooth functions the joint interaction between pairs of features and their impact on the function output can be modeled and quantified via second order partial derivatives. The Hessian of a neural network can therefore provide insight as to how interactions between pairs of features impact predictions made by the network. ReLU networks are piecewise linear and therefore have a zero Hessian almost everywhere, as a result a different approach to quantifying interactions is needed. One way to overcome this issue is to consider a smoothed approximation of the network. Prior works substitute ReLU for Softplus and compute the Hessian of the corresponding smooth network, however, the sensitivity to the $\\beta$ parameter of the Softplus activation is pronounced. Furthermore this method affords poor control over the level of smoothing. This work considers instead estimating interactions by computing the Hessian of the convolution of the network with a Gaussian - this Hessian is referred to as SmoothHess. Furthermore, using Stein's Lemma it is shown that one can efficiently estimate SmoothHess from first order information alone using Monte Carlo simulation. Non-asymptotic bounds on the error of this estimate are also provided in Theorem 3. Experimental results demonstrating improvements of this method over existing techniques, at least in the context of relatively simple problems such as MNIST and CIFAR10, are provided. ", "questions": "1. Are these techniques limited to ReLU networks or could one in fact apply them to a more general class of functions for which the first and second order partial derivatives are not defined everywhere (non-smooth) and or zero almost everywhere?\n\n2. How does one choose the covariance matrix and hence the degree of smoothing? \n\n3. Can you identify any (potentially pathological) examples in which SmoothHess really fails to model interaction effects? In short could you provide further comment on potential failure modes of the method and highlight any warnings one should be aware of when using it?\n", "limitations": "Limitations of the work, in particular the computational complexity with regard to the input dimension, are highlighted and discussed in the conclusion. A link to code is included in the supplementary material. Potential social impacts are also discussed in the supplementary material: the authors are careful to point out that as with any method used to explain decisions made by a model one needs to treat such outputs with care and due oversight, and that it is best to work in tandem with domain experts.", "rating": 7, "confidence": 3, "soundness": 3, "presentation": 3, "contribution": 3, "strengths": "This paper appears to provide a new and natural method with which to estimate the impact of interactions between pairs of features on the network output in the setting in which the network of the Hessian is zero almost everywhere. The method seems principled and outperforms existing methods on a number of benchmarks. The technical results appear correct although I confess I did not check them thoroughly. The paper is generally well written and clear. I am not an expert in this space so cannot say much to as to its broader context with respect to prior works.", "weaknesses": "One minor concern might be around how one actually chooses the covariance matrix in practice. However, choosing hyperparameters is a problem for the other methods as well it seems and furthermore this method appears to give you more direct control of the smoothing.", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "Thank you for taking the time to review our work and for the positive assessment. Below, we respond to the weaknesses and questions brought up. \n\n---\n\n## Weaknesses \n\n**W1: Choosing the covariance**\n\nAs you mention, SmoothHess gives the user more direct control over smoothing compared to the competing SoftPlus method. \n\nIn our experiment, we use isotropic covariance matrices of the form $\\sigma^2 I \\in \\mathbb{R}^{d \\times d}$. In this work we recommend a straight-forward way of choosing the level of smoothing $\\sigma$. As $d \\rightarrow \\infty$ it can be seen that the $d$-dimensional Gaussian distribution $\\mathcal{N}(0, \\sigma^2 I_{d \\times d})$ converges to a uniform distribution over the sphere of radius $\\sigma \\sqrt{d}$ [1]. in practice we have found that for the finite values of $d$ used in this paper (300,784,3072) the samples from the $d$-dimensional Gaussian $\\mathcal{N}(0, \\sigma^2 I)$ are indeed close to the sphere of radius $\\sigma \\sqrt{d}$. Thus, given one wishes to model the interactions of some network $F : \\mathbb{R}^d \\rightarrow \\mathbb{R}^c$ in a radius $\\varepsilon$-ball around a point $x_0 \\in \\mathbb{R}^d$ with SmoothHess, they may select $\\sigma = \\frac{\\varepsilon}{\\sqrt{d}}$. On the other hand, there is no clear cut way to choose the smoothing parameter $\\beta$ for SoftPlus to draw information from a specific locality. \n\nWe validate the utility of this perspective in our P-MSE experiments, shown in Table 1. In this experiment we aim to model network outputs over three neighborhoods of size $\\varepsilon = 0.25, 0.50, 1.00$ for MNIST, FMNIST and CIFAR10. SmoothHess, for which $\\sigma$ was chosen from a set of three values based on validation performance curated using the intuition above (described in lines 271-275), outperformed SoftPlus Hessian in each experiment. This is despite the fact that the best values of smoothing parameter $\\beta$ were chosen for SoftPlus Hessian from over one-hundred options, each checked on a validation set. \n\nWe additionally give an example of a use case of an alternative, non-isotropic, covariance and illustrate how such a covariance may be set. Here the goal is to focus smoothing to particular directions of input space. \n\nConsider one has a network $F: \\mathbb{R}^d \\rightarrow \\mathbb{R}^c$ and two points $x, y \\in \\mathbb{R}^d$, which they wish to quantify the interactions of $F$ between. These points may have some interesting relationship: for instance $y$ may be an adversarial example for $x$ under $F$, and one may wish to find the interactions which influence this to be the case. Let us consider the sub-network $f : \\mathbb{R}^d \\rightarrow \\mathbb{R}$, which outputs the logit (or if so desired, SoftMax probability) of the predicted class for $x$: $\\text{argmax}_{c} F(x)$. Concretely, we aim to find interactions which affect $f$ between $x$ and $y$. This can be accomplished in the following way:\n\nIntuitively, we wish to estimate SmoothHess for the point between $x$ and $y$, $z = \\frac{x + y}{2}$, and to focus smoothing along the direction $v = \\frac{y - x}{\\lVert y - x \\rVert_2}$. One may use a procedure such as Gram-Schmidt to obtain an orthonormal basis for $\\mathbb{R}^d$ for which the first vector is $v$ : $v, a_1, \\ldots, a_{d-1} \\in \\mathbb{R}^d$. Let us choose a relatively large eigenvalue of $\\sigma_l > 0$ for which to smooth along $v$ and small $\\sigma_s \\approx 0$ to smooth along the other directions. One may construct eigenvector matrix $Q = [v | a_1 | \\ldots | a_{d-1}] \\in \\mathbb{R}^{d \\times d}$ and eigenvalue matrix $\\Lambda = \\text{diag}(\\sigma_l, \\sigma_s, \\ldots, \\sigma_s) \\in \\mathbb{R}^{d \\times d}$ and the covariance may be set simply by matrix multiplication $\\Sigma = Q\\Lambda Q^T$. Finally, SmoothHess may be fit over mid-point $z$ using $\\Sigma$. One may see that this affords a higher weight to samples which are on, or close to, the line containing $x$ and $y$. \n\nConceptually this is an example of the more general approach we allude to in lines 163-165; directions of interest can be encoded with eigenvectors and their relative weighting with eigenvalues. We will include this explanation in the camera-ready Appendix.  \n\n---\n\n## Questions\n\n**Q1: Limitation to ReLU networks**\n\nThis is a great point. As we explain in global response (G1), these techniques can be applied to any function which satisfies the assumptions of Proposition 1, namely that the function is Lipschitz continuous. Thus, as long as one has access to a first-order gradient oracle, SmoothHess may be estimated for any Lipschitz continuous function. \n\n**Q2: Choosing the covariance**\n\nWe kindly direct the Reviewer to our response to W1 above. \n\n**Q3: Potential pathological examples** \n\nAs illustrated in our discussion regarding the choice of covariance, SmoothHess allows for fine-grained control over the locality of the information that is taken into account when computing interactions. Specifically, this locality is determined by the eigenvectors and eigenvalues chosen for the covariance matrix. However, in our current work, we assume that the SmoothHess smoothing is determined by a Gaussian distribution. Given one wishes to represent interactions from a specific region, they should take into account how well a Gaussian smoothing can agglomerate information from that region. Future work might investigate parameterizing the SmoothHess smoothing with different distributions, which would allow for estimating interactions over a broader range of regions for a given sample.\n\n**References** \n\n[1] Roman Vershynin. High-dimensional probability: An introduction with applications in data science, volume 47. Cambridge university press, 2018.\n"}]

Author rebuttal:
We would like to thank the reviewers for their thoughtful questions and comments, which have helped us to improve our work. 

 
---

## Recurring Questions 

Below we respond to recurring questions among the reviewers. We also give individual responses to all the reviewer questions under each reviewer’s comments

**(G1) Generality of results:** Reviewers GnHq and bZYQ raised questions about the generality of our results for non-ReLU networks or functions. As highlighted in  Proposition 1 in the paper, SmoothHess may be estimated for any Lipschitz continuous function, given one can query the gradient. We focus our presentation on ReLU networks as smoothing _is required_ in order to evaluate a non-zero Hessian. SmoothHess is still useful for Lipschitz continuous networks which have a non-zero Hessian. The reason for this is that this Hessian reflects a very small area, while SmoothHess can incorporate information from larger regions by adjusting the smoothing parameter $\Sigma$.

**(G2) Spirometry experiment:** Reviewers yoRr, FahQ, and HbQW expressed a need for clarification of the spirometry regression case-study. The purpose of the spirometry experiment is to qualitatively evaluate whether SmoothHess can be used to understand how ReLU models make predictions. FEV$_1$, a metric frequently used to evaluate lung health, is always calculated within the first 2 seconds of the spirometry curve, and _is known to be strongly affected by the presence of coughing_ [1]. In each curve, coughing may occur during, or after, the first 2 seconds. The time when coughing occurs, indicated by plateaus in the spirometry curves in Figure 3, should be an important signal for any model trained to predict FEV$_1$. 

We aim to test whether SmoothHess can detect interactions between the initial segment of the curve and segments that indicate coughing (as characterized by plateauing in the curve). To this end, we train a network to predict FEV$_1$ for spirometry samples that exhibit coughing. When we apply SmoothHess in Figure 3, we find that indeed the strongest interactions for the first 0.5 second curve segment are where plateaus occur in the curve, indicating coughing. We will add this clarification to the camera-ready revision.

---

## New Results

We also present new results in the global response pdf, which we explain below. 

**(G3) Table 1**: Reviewer HbQW asked for a comparison with the Hessian and gradient of a Swish smoothed network. Swish is a smooth activation function that may be parameterized by a smoothing value $\beta$ [2]. We compare the P-MSE of the Hessian + Gradient and Gradient of the Swish smoothed network with our previous results for the predicted class logit on MNIST, FMNIST, and CIFAR10, over three different radii $\varepsilon = 0.25, 0.50, 1.00$. Following the methodology used for SoftPlus in our work, we select the best value of $\beta$ for Swish based on the validation set performance from a set of over one-hundred options. This is done separately for Swish Hessian + Swish Gradient (Swish (H + G)) and Swish Gradient (Swish G), each dataset, and for each locality. 

Our new results indicate that the second-order Taylor expansion using SmoothHess and SmoothGrad (our method) outperforms both Swish (H + G) and Swish G in all cases. Swish (H + G) has a particularly poor performance on CIFAR10. We will include these results in our camera-ready Appendix. 

**(G4) Table 2**: Reviewers yoRr and GnHq were both interested in the ability of  SmoothHess to generalize to larger networks. We train a ResNet101 ($\approx 44.5$M parameters) on CIFAR10 and calculate P-MSE results for SmoothHess + SmoothGrad (our method), SmoothGrad, and the vanilla (unsmoothed) gradient at three radii $\varepsilon = 0.25, 0.50, 1.00$. Our method significantly outperforms both SmoothGrad and the vanilla gradient in each case. This provides evidence that SmoothHess can model interactions in large networks such as ResNet101, improving upon first-order models. 

Due to the size of ResNet101, we did not have the time to evaluate SoftPlus or Swish results for ResNet101. This is due to the expensive validation procedure (mentioned in G3) that must be used for SoftPlus and Swish. In our original P-MSE experiment, we validate over $100$ different values of SoftPlus parameter $\beta$ on a held out set, before selecting the best to use on the test set. This is because there is no clear-cut way to choose $\beta$ given a desired locality for smoothing. This is opposed to our procedure for choosing $\sigma$ for SmoothHess before computing P-MSE, which exploits well known properties of the Gaussian to curate a set of only three values of $\sigma$ for validation (as described in lines 171-177 and 269-275).

We plan to validate the best values of $\beta$ for SoftPlus and Swish on the ResNet101 and add the results to Table 2, which we will include in our camera-ready Appendix. We note that there is little reason to suspect SoftPlus or Swish will be competitive with SmoothHess. As described in (G3), SmoothHess outperforms the Swish Hessian in all cases. Further, our original results show SmoothHess outperforming SoftPlus on 17 of the 18 P-MSE experiments, with one tie. However, Reviewer FahQ inquired as to the printing precision (number of significant figures) of the reported tie, and, once checked, we found that SmoothHess achieved a lower P-MSE.

---

If any reviewers have further questions or wish us to elaborate on our responses, we would be happy to address them over this upcoming discussion period. 

---

## References

[1] Luo AZ, Whitmire E, Stout JW, Martenson D, Patel S. Automatic characterization of user errors in spirometry. Annu Int Conf IEEE Eng Med Biol Soc. 2017 Jul;2017:4239-4242. doi: 10.1109/EMBC.2017.8037792. PMID: 29060833.

[2] Ramachandran, Prajit, Barret Zoph, and Quoc V. Le. "Searching for activation functions." arXiv preprint arXiv:1710.05941 (2017).


Meta-review: Recent work has considered studying interpretability by modeling feature interactions through the Hessian of a neural network (output-to-feature). As the Hessian of ReLU networks, being piecewise-linear, (generically) vanishes, the authors propose to study ReLU networks by convolving with a Gaussian. 

Smoothing techniques of similar nature are used throughout the optimization literature. The concrete application proposed in the work, however, has been found novel and significant by the reviewers, who have also appreciated the clarity and the simplicity of the approach. Please make sure to address the important points raised by the reviewers; in particular, applicability of the proposed method beyond relatively small-scale settings.
Decision: Accept (poster)