Not Just A Black Box:
Learning Important Features Through Propagating Activation Differences

arXiv:1605.01713v3 [cs.LG] 11 Apr 2017

Avanti Shrikumar1 (avanti@stanford.edu), Peyton Greenside2 (pgreens@stanford.edu)
Anna Y. Shcherbina2 (annashch@stanford.edu), Anshul Kundaje1,3 (akundaje@stanford.edu)
1. Department of Computer Science, Stanford University, CA, USA
2. Biomedical Informatics, Stanford University, CA, USA
3. Department of Genetics, Stanford University, CA, USA

Abstract
This paper describes an older version of
DeepLIFT. See https://arxiv.org/
abs/1704.02685 for the new version. The
purported “black box” nature of neural networks
is a barrier to adoption in applications where
interpretability is essential. Here we present
DeepLIFT (Learning Important FeaTures), an
efficient and effective method for computing importance scores in a neural network. DeepLIFT
compares the activation of each neuron to its
‘reference activation’ and assigns contribution
scores according to the difference. We apply
DeepLIFT to models trained on natural images and genomic data, and show significant
advantages over gradient-based methods.

1. Introduction
As neural networks become increasingly popular, their
“black box” reputation is a barrier to adoption when interpretability is paramount. Understanding the features
that lead to a particular output builds trust with users and
can lead to novel scientific discoveries. Simonyan et al.
(2013) proposed using gradients to generate saliency maps
and showed that this is closely related to the deconvolutional nets of Zeiler & Fergus (2014). Guided backpropagation (Springenberg et al., 2014) is another variant which
only considers gradients that have positive error signal. As
shown in Figure 2, saliency maps can be substantially improved by simply multiplying the gradient with the input
signal, which corresponds to a first-order Taylor approximation of how the output would change if the input were
set to zero; as we show, the layer-wise relevance propaProceedings of the 33 rd International Conference on Machine
Learning, New York, NY, USA, 2016. JMLR: W&CP volume
48. Copyright 2016 by the author(s).

gation rules described in Bach et al. (2015) reduce to this
approach, assuming bias terms are included in the denominators.
Gradient-based approaches are problematic because activation functions such as Rectified Linear Units (ReLUs) have
a gradient of zero when they are not firing, and yet a ReLU
that does not fire can still carry information (Figure 1).
Similarly, sigmoid or tanh activations are popular choices
for the activation functions of gates in memory units of recurrent neural networks such as GRUs and LSTMs (Chung
et al., 2014; Hochreiter & Schmidhuber, 1997), but these
activations have a near-zero gradient at high or low inputs
even though such inputs can be very significant.

Figure 1. Simple network with inputs x1 and x2 that have reference values of 0. When x1 = x2 = −1, output is 0.1 but the
gradients w.r.t x1 and x2 are 0 due to inactive ReLU y (which has
activation of 2 under reference input). By comparing activations
to their reference values, DeepLIFT assigns contributions
to the

output of (0.1 − 0.5) 13 to x1 and (0.1 − 0.5) 23 to x2 .

We present DeepLIFT, a method for assigning feature importance that compares a neuron’s activation to its ‘reference’, where the reference is the activation that the neuron
has when the network is provided a ‘reference input’ (the
reference input is defined according to what is appropriate for the task at hand). This addresses the limitation of
gradient-based approaches because the difference from the
reference may be non-zero even when the gradient is zero.

2. DeepLIFT Method
We denote the contribution of neuron x to neuron y as Cxy .
Let the activation of a neuron n be denoted as An . Further,
let the reference activation of neuron n be denoted A0n , and
let the An − A0n be denoted as δn . We define our contribu-

This paper describes an older version of DeepLIFT; see https://arxiv.org/abs/1704.02685 for the new version

tions Cxy to satisfy the following properties.

low:

2.1. Summation to δ

X Cxy
Cyt
δy

Cxt =

For any set of neurons S whose activations are minimally
sufficient to compute the activation of y (that is, if we know
the activations of S, we can compute the activation of y,
and there is no set S 0 ⊂ S such that S 0 is sufficient to
compute the activation of y - in layman’s terms, S is a full
set of non-redundant inputs to y), the following property
holds:
X
Csy = δy
(1)

y∈Ox

X
X Cxy
Cxy myt
(myt δy ) =
δy

mxt δx =

(5)

y∈Ox

y∈Ox

X Cxy
X
mxy myt
myt =
δx

mxt =

y∈Ox

y∈Ox

In the equations below, Iy denotes the set of inputs of y.

s∈S

2.3.1. A FFINE FUNCTIONS
That is, the sum over all the contributions of neurons in S
to y equals the difference-from-reference of y.

Let


Ay = 

2.2. Linear composition

X Cxo
Coy
δo

(2)

o∈Ox

In layman’s terms, each neuron ‘inherits’ a contribution
through its outputs in proportion to how much that neuron
contributes to the difference-from-reference of the output.

(6)

Then mxy = wxy
Proof. We show that δy =

x∈Iy

mxy δx .

x∈Iy


=

(7)


X

wxy A0x  + b +

x∈Iy

We show that the contributions as defined above can be
computed using the following rules (which can be implemented to run on a GPU). The computation is reminiscent
of the chain rule used during gradient backpropagation,
as equation 2 makes it possible to start with contribution
scores of later layers and use them to find the contribution
scores of preceding layers. To avoid issues of numerical
stability when δn for a particular neuron is small, rather
than computing the contribution scores explicitly, we instead compute multipliers mxy that, when multiplied with
the difference-from-reference, give the contribution:
mxy δx = Cxy

(3)

Let t represent the target neuron that we intend to compute
contributions to, and let Ox represent the set of outputs of
x. We show that:
X

P

Using the fact that An = A0n + δn , we have:


X
(A0y + δy ) = 
wxy (A0x + δx ) + b

2.3. Backpropagation Rules

mxt =

wxy Ax  + b

x∈Iy

Let Ox represent the output neurons of x. The following
property holds:
Cxy =


X

mxy myt

(4)

y∈Ox

The equation above follows from the linear composition
property and the definition of the multipliers, as proved be-

X

wxy δx

x∈Iy

We also note that the reference activation A0y can be found
as follows:


X
A0y = 
wxy A0x  + b
(8)
x∈Iy

Thus, canceling out
δy =

A0y

X

yields:
X

wxy δx =

x∈Iy

mxy δx

(9)

x∈Iy

2.3.2. M AX OPERATION
We consider the case of max operation such as a maxpool:
Ay = max Ax

(10)

x∈Iy

Then we have:
mxy = 1{Ax = Ay }

δy
δx

(11)

Where 1{} is the indicator function. If a symbolic computation package is used, then the gradient of y with respect
to x can be used in place of 1{Ax = Ay }.

This paper describes an older version of DeepLIFT; see https://arxiv.org/abs/1704.02685 for the new version

2.4. A note on final activation layers

Proof.
X

X

mxy δx =

x∈y

1{Ax = Ay }

x∈y

=

X

δy
δx

!
δx
(12)

1{Ax = Ay }δy = δy

x∈y

2.3.3. M AXOUT UNITS
A maxout function has the form
!
X
n
i
Ay = max
wxy Ax + bi
i=1

(13)

x

i.e. it is the max over n affine functions of the input vector ~x. For a given vector of activations A~x of the inputs,
we split A~x − A~0x into segments such that over each segment s, a unique affine function dominates the maxout and
the coefficient of an individual input x over that segment is
w(s)xy . Let l(s) denote the fraction of A~x −A~0x in segment
s. We have:
X
mxy =
l(s)w(s)xy
(14)
s

Intuitively speaking, we simply split the piecewise-linear
maxout function into regions where it is linear, and do a
weighted sum of the coefficients of x in each region according to how much of A~x − A~0x falls in that region.
2.3.4. OTHER ACTIVATIONS
The following choice for mxy , which is the same for all
inputs to y, satisfies summation-to-delta:
mxy = P

δy

x0 ∈I

y

(15)

δ x0

This rule may be used for nonlinearities like ReLUs, PReLUs, sigmoid and tanh (where y has only one input). Situations where the denominator is near zero can be handled
by applying L’hopital’s rule, because by definition:
X
δy → 0 as
δx → 0
(16)
x∈Iy

2.3.5. E LEMENT- WISE PRODUCTS

Activation functions such as a softmax or a sigmoid have a
maximum δ of 1.0. Due to the summation to δ property, the
contribution scores for individual features are lower when
there are several redundant features present. As an example, consider At = σ(Ay ) (where sigma is the sigmoid
transformation) and Ay = Ax1 + Ax2 . Let the default activations of the inputs be A0x1 = A0x2 = 0. When x1 = 100
and x2 = 0, we have Cx1 t = 0.5. However, when both
x1 = 100 and x2 = 100, we have Cx1 t = Cx2 t = 0.25.
To avoid this attenuation of contribution in the presence of
redundant inputs, we can use the contributions to y rather
than t; in both cases, Cx1 y = 100.
2.5. A note on Softmax activation
Let t1 , t2 ...tn represent the output of a softmax transformation on the nodes y1 , y2 ...yn , such that:
eAyi
Ati = Pn
Ay0
i
i0 =1 e

Here, Ay1 ...Ayn are affine functions of their inputs. Let
x represent a neuron that is an input to Ay1 ...Ayn , and
let wxyi represent the coefficient of Ax in Ayi . Because
Ay1 ...Ayn are followed by a softmax transformation, if
wxyi is the same for all yi (that is, x contributes equally to
all yi ), then x effectively has zero contribution to Ati . This
can be observed by substituting Ayi = wxyi Ax + ryi in the
expression for Ati and canceling out ewxyi Ax (here, ryi is
the sum of all the remaining terms in the affine expression
for Ayi )
ewxyi Ax +ryi
eAyi
=
Ati = Pn
P
n
Ay0
wxy 0 Ax +ry 0
i
i
i
i0 =1 e
i0 =1 e
= Pn

ewxyi Ax +ryi

i0 =1

ewxyi Ax +ryi0

Ay =

+ δy =

(A0x1

+

δx1 )(A0x2

+ δ x2 )

(17)

We have:
δy = (A0x1 + δx1 )(A0x2 + δx2 ) − (A0x1 A0x2 )
= A0x1 δx2 + A0x2 δx1 + δx1 δx2




δx2
δ x1
0
0
= δ x 1 Ax 2 +
+ δ x 2 Ax 1 +
2
2
Thus, viable choices for the multipliers are mx1 y =
0.5δx2 and mx2 y = A0x1 + 0.5δx1

(18)

eryi
= Pn
ry 0
i
i0 =1 e

(20)

As mentioned in the previous subsection, in order to avoid
attenuation of signal for highly confident predictions, we
should compute Cxyi rather than Cxti . One way to ensure
that Cxyi is zero if wxyi is the same for all yi is to meannormalized the weights as follows:

Consider the function:
A0y

(19)

w̄xyi = wxyi −

n
1X
wxyi0
n 0

(21)

i =1

This transformation will not affect the output of the softmax, but will ensure that the DeepLIFT scores are zero
when a particular node contributes equally to all softmax
classes.
2.6. Weight normalization for constrained inputs

A0x2

+

Let y be a neuron with
P some subset of inputs Sy that are
constrained such that x∈Sy Ax = c (for example, one-hot

This paper describes an older version of DeepLIFT; see https://arxiv.org/abs/1704.02685 for the new version

P
encoded input satisfies the constraint x∈Sy Ax = 1, and
a convolutional neuron operating on one-hot encoded rows
has one constraint per column that it sees). Let the weights
from x to y be denoted wxy and let by be the bias of y. It is
advisable to use normalized weights w̄xy = wxy − µ and
bias b̄y = by + cµ, where µ is the mean over all wxy . We
note that this maintains the output of the neural net because,
for any constant µ:
X


Ax (w̄xy − µ) + (by + cµ)
X
 X

=
Ax wxy −
Ax µ + (by + cµ)
X

=
Ax wxy − cµ + (by + cµ)
X

=
Ax wxy + by

Ay =

(22)

The normalization is desirable because, for affine functions, the multipliers mxy are equal to the weights wxy
and are thus sensitive to µ. To take the example of a convolutional neuron operating on one-hot encoded rows: by
mean-normalizing wxy for each column in the filter, one
can ensure that the contributions Cxy from some columns
are not systematically overestimated or underestimated relative to the contributions from other columns.

3. Results
3.1. Tiny ImageNet
A model with the VGG16 (Long et al., 2015) architecture
was trained using the Keras framework (Chollet, 2015) on a
scaled-down version of the Imagenet dataset, dubbed ‘Tiny
Imagenet’. The images were 64 × 64 in dimension and
belonged to one of 200 output classes. Results shown in
Figure 2; the reference input was an input of all zeros after
preprocessing.
3.2. Genomics
We apply DeepLIFT to models trained on genomic sequence. The positive class requires that the DNA patterns ’GATA’ and ’CAGATG’ appear in the length-200 sequence together. The negative class has only one of the
two patterns appearing once or twice. Outside the core patterns (which were sampled from a generative model) we
randomly sample the four bases A, C, G and T. A CNN
was trained using the Keras framework (Chollet, 2015) on
one-hot encoded sequences with 20 convolutional filters of
length 15 and stride 1 and a max pool layer of width and
stride 50, followed by two fully connected layers of size
200. PReLU nonlinearities were used for the hidden layers. This model performs well with auROC of 0.907. The
misclassified examples primarily occur when one of the
patterns erroneously arises in the randomly sampled background.

Figure 2. Comparison of methods. Importance scores for RGB
channels were summed to get per-pixel importance. Left-toright: original image, absolute value of the gradient (similar to
Simonyan et al. which used the two-norm across RGB rather than
the sum, and which is related to both Zeiler et al. and Springenberg et al.), positive gradient*input (Taylor approximation, equivalent to Layer-wise Relevance Propagation in Samek et al. but
masking negative contributions), and positive DeepLIFT.

We then run DeepLIFT to assign an importance score to
each base in the correctly predicted sequences. The reference input is an input of all zeros post weight-normalization
(see 2.6) of the first convolutional layer (after weight normalization, the linear activation of a convolutional neuron
for an input of all zeros is the bias, which is the same as the
average activation across all four bases at each position).
We compared the results to the gradient*input (Figure 3).

Figure 3. DeepLIFT scores (top) and gradient*input (bottom) are
plotted for each position in the DNA sequence and colored by
the DNA base (due to one-hot encoding, input is either 1 or 0;
gradient*input is equivalent to taking the gradient for the letter
that is actually present). DeepLIFT discovers both patterns and
assigns them large importance scores. Gradient-based methods
miss the GATA pattern.

4. Discussion
Prevailing feature importance methods such as the saliency
maps of Simonyan et al., the deconvolutional nets of Zeiler
et al. and the guided backpropagation of Springenberg et
al. are variants of computing gradients. As shown in Figure 1, this can give misleading results when the local gradient is zero. DeepLIFT instead considers the deviation
from a neuron’s reference activity. This makes it capable
of handling RNN memory units gated by activations that
have vanishing gradients (eg: sigmoid, tanh).
Layer-wise Relevance Propagation (LRP), proposed by

This paper describes an older version of DeepLIFT; see https://arxiv.org/abs/1704.02685 for the new version

Bach et al., does not obviously rely on gradients; however,
as we show, if all activations are piecewise linear, LRP reduces to gradient*input (a first-order Taylor approximation
of the change in output if the input is set to zero). If all reference activations are zero (as happens when all bias terms
are zero and all reference input values are zero), DeepLIFT
and LRP give similar results (except that by computing
contributions using multipliers, DeepLIFT circumvents the
numerical stability problems that LRP faces). In practice,
biases are often non-zero, which is why DeepLIFT produces superior results (Figures 2 & 3).

Let us now consider what happens when there are two filtering operations applied sequentially. Let Rik denote the
relevance inherited by neuron i in layer l from neuron k in
layer l + 2, passing through the neurons in layer l + 1. We
have:

4.1. Equivalence of gradient*input to Layer-wise
Relevance Propagation

Thus, we see that denominator aj

Unpooling: “The backwards signal is redirected proportionally onto the location for which the activation was
recorded in the forward pass”: This is trivially the same
as gradient*input, because the gradient*input will be zero
for all locations which do not activation the pooling layer,
and equal to the output for the location that does.
Filtering: We consider the first rule described in Samek et
(l) (l,l+1)
al., where zij = ai wij
is the weighted activation of
neuron i onto neuron j in the next layer, and l is the index
of the layer:
X
zij
(l+1)
(l)
Rj
(23)
Ri =
P0
P
0 j +  sign(
0j )
z
z
0
i
i
i
i
j
The term involving  P
is included to avoid issues of numer0
ical instability when i zi0 j is near zero. The second rule
described in Samek et al. is another variant designed to address the problem of numerical instability. We show that
gradient*input gives the exact result as  → 0 (i.e. it solves
the issue of numerical instability altogether).
(l,l+1)

Dropping the term for  and substituting zij = ai wij
we have:
(l)

(
Ri l)

=

X
j

,

(l,l+1)

ai wij
(l+1)
P0 (l) (l,l+1) Rj
i ai0 wi0 j

(l,l+1) (l+1) (l+1,l+2)
X a(l)
aj
wjk
i wij
(l+1)

j

=

(l+1)

aj

ak

(l,l+1) (l+1,l+2)
X a(l)
wjk
i wij
(l+1)

ak

j

(24)

Assuming the bias term is included (which would be necessary for the conservation property described in Bach et al.
to hold), the denominator is simply the activation of neuron
j, i.e.:
(l,l+1)
X a(l)
i wij
(
(l+1)
Ri l) =
Rj
(25)
(l+1)
a
j
j

(l+2)

Rk

(26)
(l+2)
Rk

(l+1)

We show when all activations are piecewise linear and bias
terms are included in the calculation, the Layer-wise Relevance Propagation (LRP) of Bach et al., reduces to gradient*input. We refer to Samek et al. (2015) for the concise
description of LRP:

(l)

(l)

Rik =

for the intermediate

(l) (l,l+1) (l+1,l+2)
layer cancelled out, leaving us with ai wij
wjk
,
(l,l+1) (l+1,l+2)
(l+1)
where wij
wjk
is the gradient of ak
with re(l)
spect to ai . The only term left in the denominator is the
(l+1)
activation of the last layer, ak
; if we set the relevance of

neurons in the final layer to be equal to their own activation,
(l+2)
then Rk
(assuming k is the last layer) would cancel out
(l+1)
ak
in the denominator, leaving us with:
X (l) (l,l+1) (l+1,l+2)
(l)
Rik =
ai wij
wjk
(27)
j
(l)

Which is simply equal to the activation ai multiplied by
(l)
the gradient of ak with respect to ai . In situations where
the relevance of the last layer is not the same as its activation (which may happen if there is a nonlinear transformation such as a sigmoid, as a sigmoid output of 0.5 occurs when the input is 0), one can simply compute gradient*input with respect to the linear term before the final
nonlinearity (which is what we did; for softmax layers, we
apply the normalization described in 2.5).
Nonlinearity: “The backward signal is simply propagated
onto the lower layer, ignoring the rectification operation”:
While this is not obviously the same as gradient*input, it
should be noted that when a rectified linear unit is inactive,
it has an activation of zero and the rule for filtering (described above) would assign it zero importance. Furthermore, when the rectified linear unit is active, its gradient is
1. Thus, when the unit is inactive, gradient*input is 0 and
LRP assigns 0 signal; when a unit is active, gradient*input
is equal to the output and LRP assigns all signal. The two
approaches converge.

5. Author contributions
AS & PG conceived of DeepLIFT. AS implemented
DeepLIFT in software. PG led application to genomics.
AYS led application to Tiny Imagenet. AK provided guidance and feedback. AS, PG, AYS & AK prepared the
manuscript.

This paper describes an older version of DeepLIFT; see https://arxiv.org/abs/1704.02685 for the new version

References
Bach, Sebastian, Binder, Alexander, Montavon, Grégoire,
Klauschen, Frederick, Müller, Klaus-Robert, and
Samek, Wojciech. On Pixel-Wise explanations for
Non-Linear classifier decisions by Layer-Wise relevance
propagation. PLoS One, 10(7):e0130140, 10 July 2015.
Chollet, Franois. Keras, 2015.
Chung, Junyoung, Gulcehre, Caglar, Cho, Kyunghyun, and
Bengio, Yoshua. Empirical evaluation of gated recurrent
neural networks on sequence modeling. 11 December
2014.
Hochreiter, S and Schmidhuber, J. Long short-term memory. Neural Comput., 9(8):1735–1780, 15 November
1997.
Simonyan, Karen, Vedaldi, Andrea, and Zisserman, Andrew. Deep inside convolutional networks: Visualising
image classification models and saliency maps. 20 December 2013.
Springenberg, Jost Tobias, Dosovitskiy, Alexey, Brox,
Thomas, and Riedmiller, Martin. Striving for simplicity: The all convolutional net. 21 December 2014.
Zeiler, Matthew D and Fergus, Rob. Visualizing and understanding convolutional networks. In Computer vision–
ECCV 2014, pp. 818–833. Springer, 2014.

