Improving performance of deep learning models with axiomatic attribution priors and expected gradients

TL;DR

This paper introduces a new attribution method called expected gradients that satisfies key interpretability axioms, enabling the use of attribution priors to improve deep learning model performance across various real-world datasets.

Contribution

The paper presents a novel attribution method, expected gradients, that satisfies interpretability axioms and enhances attribution priors for better model performance.

Findings

Expected gradients satisfy key interpretability axioms.

Combining attribution priors with expected gradients improves performance.

Results are consistent across image, gene expression, and healthcare datasets.

Abstract

Recent research has demonstrated that feature attribution methods for deep networks can themselves be incorporated into training; these attribution priors optimize for a model whose attributions have certain desirable properties -- most frequently, that particular features are important or unimportant. These attribution priors are often based on attribution methods that are not guaranteed to satisfy desirable interpretability axioms, such as completeness and implementation invariance. Here, we introduce attribution priors to optimize for higher-level properties of explanations, such as smoothness and sparsity, enabled by a fast new attribution method formulation called expected gradients that satisfies many important interpretability axioms. This improves model performance on many real-world tasks where previous attribution priors fail. Our experiments show that the gains from combining higher-level attribution priors with expected gradients attributions are consistent across image, gene expression, and health care data sets. We believe this work motivates and provides the necessary tools to support the widespread adoption of axiomatic attribution priors in many areas of applied machine learning. The implementations and our results have been made freely available to academic communities.

Figures

Tables

Table 1. Table 1: Results from benchmark software on synthetic data with correlated features. Larger numbers mean a better feature attribution method for all metrics. Metric abbreviations are: K (Keep), R (Remove); P (Positive), N (Negative), A (Absolute); M (Mean masking), R (Resample masking), I (Impute masking). For example, KPM corresponds to the “Keep Positive with Mean masking” metric. Each model is trained on 900 900 900 samples and tested using 100 100 100 samples. See Supplement Section LABEL:sec:supp_benchmark for details.

Method	KPM	KPR	KPI	KNM	KNR	KNI	KAM	KAR	KAI
Expected Grad.	3.731	3.800	3.973	3.615	3.551	3.873	0.906	0.903	0.919
Integrated Grad.	3.667	3.736	3.920	3.543	3.476	3.808	0.905	0.899	0.920
Gradients	0.096	0.122	0.099	0.076	-0.112	0.052	0.838	0.823	0.887
Random	0.033	0.106	0.077	-0.012	-0.093	-0.053	0.593	0.583	0.715
Method	RPM	RPR	RPI	RNM	RNR	RNI	RAM	RAR	RAI
Expected Grad.	3.612	3.575	3.525	3.759	3.830	3.683	0.897	0.885	0.880
Integrated Grad.	3.539	3.503	3.365	3.687	3.754	3.543	0.872	0.859	0.822
Gradients	0.035	-0.098	-0.020	0.110	0.105	0.108	0.729	0.712	0.616
Random	-0.053	-0.100	-0.106	0.034	0.092	0.111	0.400	0.400	0.275

Equations

$$\theta ={\text{argmin}}_{\theta }\mathcal{L}(\theta ;X,y)+{\lambda }^{\prime }{\Omega }^{\prime }(\theta ).$$

$$\theta ={\text{argmin}}_{\theta }\mathcal{L}(\theta ;X,y)+\lambda \Omega (\Phi (\theta ,X)),$$

$${\Omega }_{\text{pixel}}(\Phi (\theta ,X))=\sum _{\ell }\sum _{i,j}|{\varphi }_{i+1,j}^{\ell }-{\varphi }_{i,j}^{\ell }|+|{\varphi }_{i,j+1}^{\ell }-{\varphi }_{i,j}^{\ell }|.$$

$${\Omega }_{\text{graph}}(\Phi (\theta ,X))=\sum _{i,j}{W}_{i,j}({\bar{\varphi }}_i-{\bar{\varphi }}_j{)}^2={\bar{\varphi }}^T{L}_G\bar{\varphi },$$

$${\Omega }_{\text{sparse}}(\Phi (\theta ,X))=-\frac{{\sum }_{i=1}^p{\sum }_{j=1}^p|{\bar{\varphi }}_i-{\bar{\varphi }}_j|}{n{\sum }_{i=1}^p{\bar{\varphi }}_i}=-2G(\bar{\varphi }),$$

$${\text{ExpectedGradients}}_i(x)$$

$$\displaystyle=\int_{x^{\prime}}\Bigg{(}(x_{i}-x_{i}^{\prime})\times\int_{\alpha=0}^{1}\frac{\delta f(x^{\prime}+\alpha(x-x^{\prime}))}{\delta x_{i}}d\alpha\Bigg{)}p_{D}(x^{\prime})dx^{\prime},$$

$$\displaystyle=\mathop{\mathbb{E}}_{x^{\prime}\sim D,\alpha\sim U(0,1)}\bigg{[}(x_{i}-x_{i}^{\prime})\times\frac{\delta f(x^{\prime}+\alpha\times(x-x^{\prime}))}{\delta x_{i}}\bigg{]}.$$

$$\theta ={\text{argmin}}_{\theta }\mathcal{L}(\theta ;X,y)+{\lambda }^{\prime \prime }||A\odot \frac{\partial \mathcal{L}}{\partial X}|{|}_F^2.$$

$$\theta ={\text{argmin}}_{\theta }\mathcal{L}(\theta ;X,y)+\lambda \Omega (\Phi (\theta ,X)),$$

$${\text{IntegratedGradients}}_i(x,x^{\prime }):=(x_i-x_i^{\prime })\times {\int }_{\alpha =0}^1\frac{\delta f(x^{\prime }+\alpha (x-x^{\prime }))}{\delta x_i}d\alpha ,$$

$$\textrm{ExpectedGradients}_{i}(x)\vcentcolon=\int_{x^{\prime}}\Bigg{(}(x_{i}-x_{i}^{\prime})\times\int_{\alpha=0}^{1}\frac{\delta f(x^{\prime}+\alpha(x-x^{\prime}))}{\delta x_{i}}d\alpha\Bigg{)}p_{D}(x^{\prime})dx^{\prime},$$

$$\textrm{ExpectedGradients}_{i}(x)\vcentcolon=\\ \mathop{\mathbb{E}}_{x^{\prime}\sim D,\alpha\sim U(0,1)}\bigg{[}(x_{i}-x_{i}^{\prime})\times\frac{\delta f(x^{\prime}+\alpha\times(x-x^{\prime}))}{\delta x_{i}}\bigg{]}.$$

$${\text{argmin}}_{\theta }||f_{\theta }(X)-y|{|}_2^2={\text{argmax}}_{\theta }\text{exp}(-||f_{\theta }(X)-y|{|}_2^2)={\theta }_{MLE},$$

$${\text{argmin}}_{\theta }||f_{\theta }(X)-y|{|}_2^2+\lambda ||\theta |{|}_2^2={\text{argmax}}_{\theta }\text{exp}(-||f_{\theta }(X)-y|{|}_2^2)\text{exp}(-\lambda ||\theta |{|}_2^2)={\theta }_{MAP},$$

$${\Omega }_{\text{pixel}}(\Phi (\theta ,X))=\sum _{\ell }\sum _{i,j}|{\varphi }_{i+1,j}^{\ell }-{\varphi }_{i,j}^{\ell }|+|{\varphi }_{i,j+1}^{\ell }-{\varphi }_{i,j}^{\ell }|,$$

$${\Omega }_{\text{graph}}(\Phi (\theta ,X))=\sum _{i,j}{W}_{i,j}({\bar{\varphi }}_i-{\bar{\varphi }}_j{)}^2={\bar{\varphi }}^T{L}_G\bar{\varphi }.$$

$${\Omega }_{\text{sparse}}(\Phi (\theta ,X))=-\frac{{\sum }_{i=1}^p{\sum }_{j=1}^p|{\bar{\varphi }}_i-{\bar{\varphi }}_j|}{n{\sum }_{i=1}^p{\bar{\varphi }}_i}=-2G(\bar{\varphi }),$$

$$\prod _{i=1}^p\prod _{j=1}^p\text{exp}\left(\frac{1}{{\sum }_{i=1}^p{\bar{\varphi }}_i}|{\bar{\varphi }}_i-{\bar{\varphi }}_j|\right),$$

$$\mathcal{L}(w;X,y)=\Vert Xw-y{\Vert }_2^2+{\lambda }^{\prime }\Vert w{\Vert }_1+{\nu }^{\prime }{w}^T{L}_{G}w,$$

$$\mathcal{L}(\theta ;X,y)=\Vert f_{\theta }(X)-y{\Vert }_2^2+\lambda \Vert \theta {\Vert }_1$$

$$\mathcal{L}(\theta ;X)={\Omega }_{graph}(\Phi (\theta ,X))=\nu {\bar{\varphi }}^T{L}_G\bar{\varphi }$$

Full text

One sentence summary: We introduce a method for using axiomatic feature attributions to encourage a wide variety of desirable behaviors in deep models.

1 Introduction

Recent work on interpreting machine learning (ML) models focuses on feature attribution methods. Given an input datum, a model, and a prediction, such methods assign a number to each input feature that represents how important the feature was for making the prediction. Current research is also investigating the axioms that attribution methods should satisfy [23, 47, 45, 7] and how they provide insight into model behavior [25, 24, 36, 54]. Feature attribution methods often reveal problems in a model or dataset. For example, a model may place too much importance on undesirable features, rely on many features when sparsity is desired, or be sensitive to high frequency noise. In such cases, humans often have a prior belief about how a model should treat input features but find it difficult to mathematically encode this prior for neural networks in terms of the model parameters.

Ross et al. [34] introduced the first instance of what we call an attribution prior: if it is possible for explanations to reveal problems in a model, constraining the model’s explanations during training can help the model avoid such problems. Given a binary indicator of whether each feature should or should not be important for making predictions on each sample in the dataset, their method penalizes the gradients of unimportant features. However, two drawbacks limit this method’s applicability to real-world problems. First, gradients do not satisfy the same theoretical guarantees as modern feature attribution methods. This leads to well-known problems such as saturation: operations, like ReLUs and sigmoids, which have large flat “saturated” regions, can lead to 0 gradient attribution even for important features [47]. Second, it can be difficult to specify which features should be important in a binary manner. A very recent publication successfully uses Ross et al’s gradient-based prior as part of a human-in-the-loop strategy to improve model generalization performance and user trust, as well as contributing their own model-agnostic method for penalizing feature importances [38]. Such results create a clear synergy with our study, which improves the quality of calculated importances and develops new forms of attribution priors. This has the potential to greatly expand both the number of ways that a human-in-the-loop can influence deep models and the precision with which they can do so.

More recent work discusses the need for priors that incorporate human intuition in order to develop robust and interpretable models [13]. Still, it remains challenging to encode priors such as “have smoother attributions across an image” or “treat this group of features similarly” by penalizing a model’s input gradients or parameters. Some recent attribution priors have proposed regularizing integrated gradients attributions [19, 5]. While promising, this work suffers from three major weaknesses: it does not clearly demonstrate improvements over gradient-based attribution priors, it penalizes attribution deviation from a target value rather than encoding sophisticated priors such as those we mention above, and it imposes a large computational cost by training with tens to hundreds of reference samples per batch.

A contemporary method called contextual decomposition explanation penalization (CDEP) uses a framework similar to attribution priors and penalizes explanations generated by the contextual decomposition (CD) method [31]. CDEP claims to greatly increase computational efficiency compared to the method presented in this paper; however, we demonstrate here that these observed performance differences were actually the result of errors in the re-implementation of the method (Supplement Section LABEL:sec:supp_cdep). In many cases, our method actually trains faster than CDEP. Our approach also maintains several advantages, including a much more flexible attribution method that works with a wider array of models (any model with gradients) and penalties (any differentiable function of the attributions for each feature). For further discussion of related work, including a discussion of cases where CDEP should work well, see Supplement Section LABEL:sec:_supp_related_work.

The main contribution of this work is a broadened interpretation of attribution priors that includes any case in which the training objective incorporates differentiable functions of a model’s axiomatic feature attributions. This can be seen as a generalization of gradient-based regularization [18, 34, 53, 14, 35] and it can be used to encode meaningful domain knowledge more effectively than existing methods. Whereas previous attribution priors generally took the form of “encourage feature $i$’s attribution to be near a pre-determined target value,” the priors we present here consider relative importance among multiple features and do not require pre-determined target values for any feature’s attribution. Specifically, we introduce an image prior enforcing that neighboring pixels have similar attributions, a graph prior for biological data enforcing that related genes have similar attributions, and a sparsity prior enforcing that a few features have large attributions while all others have near-zero ones.

We also introduce a new general-purpose feature attribution method to enforce these priors, expected gradients, which eliminates a hyperparameter choice required by integrated gradients [47]. Further, it can be efficiently regularized under an attribution prior due to its formulation as an expectation, which naturally lends itself to batched estimates of the attribution. Since these attributions are used not only to interpret trained models, but also as part of the training objective itself, it is essential to guarantee that the attributions will be of high quality. We therefore show that our attribution method satisfies important interpretability axioms.

Across three different prediction tasks, we show that training with expected gradients outperforms training with previous, more limited versions of attribution priors. On images, our image prior produces a model that is more interpretable and generalizes better to noisy data. On gene expression data, our graph prior reduces prediction error and better captures biological signal. Finally, on a patient mortality prediction task, our sparsity prior yields a sparser model and improves performance when learning from limited training data.

2 Results

2.1 Attribution priors are a flexible framework for encoding domain knowledge.

Let $X\in\mathbb{R}^{n\times p}$ denote a dataset with labels $y\in\mathbb{R}^{n\times o}$, where $n$ is the number of samples, $p$ is the number of features, and $o$ is the number of outputs. In standard deep learning, we find optimal parameters $\theta$ by minimizing loss, with a regularization term $\Omega^{\prime}(\theta)$ weighted by $\lambda^{\prime}$ on the parameters:

[TABLE]

Attribution priors involve a model’s attributions, represented by the matrix $\Phi(\theta,X)$, where each entry $\phi_{i}^{\ell}$ is the importance of feature $i$ in the model’s output for sample $\ell$. The attribution prior is a scalar-valued penalty function of the feature attributions $\Omega(\Phi(\theta,X))$, which represents a log-transformed prior probability distribution over possible attributions ($\lambda$ is the regularization strength). The attribution prior is modular and agnostic to the particular attribution method. This results in the optimization:

[TABLE]

where the standard regularization term has simply been replaced with an arbitrary, differentiable penalty function on the feature attributions.

While feature attributions have previously been used in training (more details in Methods 1) [34, 19], our approach offers two novel components: First, we demonstrate that calculating $\Phi$ with attribution methods that satisfy previously-established interpretability axioms improves performance (see Section 2.2 and Methods 2 for further discussion of interpretability axioms). Second, rather than simply encouraging each feature’s attribution to be near a target value as in previous work, we enforce high-level priors over the relationships between features. In image data, we use a Laplace 0-mean prior on the difference between attributions of adjacent pixels, which encourages a low total variation (high smoothness) of attributions:

[TABLE]

In gene expression data, we use a Gaussian 0-mean prior on the difference between mean absolute attributions $\bar{\phi}_{i}$ of functionally related genes, which encourages such similar genes to have similar attributions:

[TABLE]

where $W_{i,j}$ is the weight of connection between two genes in a biological graph, and $L_{G}$ is the graph Laplacian. Finally, in health data where sparsity is desired, we use a prior on the Gini coefficient of the mean absolute attributions $\bar{\phi}_{i}$, which encourages a small number of features to have a large percentage of the total attribution while others are near-zero:

[TABLE]

where $G$ is the Gini coefficient. None of these priors require specifying target values for features, and all improve performance over simpler baselines. For more details on our priors see Methods 3, and for more details on previous attribution priors, see Methods 1.

2.2 Expected gradients outperforms other attribution methods.

Attribution priors involve using feature attributions not just as a post-hoc analysis method, but as a key part of the training objective. Thus, it is essential to guarantee as much as possible that the attribution method used will produce high-quality attributions. We propose an axiomatic feature attribution method called expected gradients, which avoids problems with existing methods and is naturally suited to being incorporated into training. Expected gradients extends the integrated gradients method [47], and like integrated gradients, satisfies a variety of desirable interpretability axioms such as completeness (the feature attributions sum to the output for a given sample) and implementation invariance (the attributions are identical for any of the infinite possible implementations of the same function). Because these methods satisfy completeness, they are not subject to problems with input saturation, while methods such as gradients that do not satisfy completeness will suffer from input saturation. Because these methods satisfy implementation invariance, they are straight-forward to practically apply to any differentiable model, regardless of specific network architectures (see Methods 2 for an extended discussion of the interpretability axioms satisfied by expected gradients).

Integrated gradients generates feature attributions by explaining the difference between the model’s output for a sample and its output for a reference sample. Choosing this reference is difficult. For example, in image tasks, the image of all zeros is often chosen as a baseline, but doing so implies that black pixels will not be highlighted as important. This problem can be solved by integrating gradients over multiple references. However, this requires multiple integrals, is expensive in terms of time and memory, and may be prohibitive to calculate for every sample in every epoch of training. Our method collapses these multiple integrals using sampling to yield accurate attributions with multiple references that can be calculated quickly (here, $x$ is the sample, $x^{\prime}$ is a reference, and $D$ is the data distribution):

[TABLE]

If the attribution function $\Phi$ in our attribution prior $\Omega(\Phi(\theta,X))$ is integrated gradients, regularizing $\Phi$ would require hundreds of extra gradient calls every training step (the original integrated gradients paper [47] recommends 20 to 300 gradient calls to compute attributions). This makes training with integrated gradients prohibitively slow – in fact, [19] find that using integrated gradients can take up to 30 times longer than standard training even when only back-propagating gradients through part of the network. Even more gradient calls would be necessary if multiple references were used, as proposed above. However, most deep learning models today are trained using some variant of batch gradient descent, where the gradient of a loss function is approximated over many training steps using mini-batches of data. We can use a batch training procedure to approximate expected gradients as well. During training, we let $k$ be the number of samples we draw to compute expected gradients for each mini-batch. Remarkably, as small as $k=1$ suffices to regularize the explanations because of the averaging effect of the expectation formulation over many training samples. This choice of $k$ uses every sample in the training set as a reference over the course of an epoch, with only one additional gradient call per training step. This lets us efficiently regularize expected gradients using the entire training dataset as the reference – as opposed to the single reference in integrated gradients – with far fewer gradient calls per epoch (see Methods 2 for more detail).

We compare expected gradients to other feature attribution methods for 18 benchmarks on synthetic data introduced in [25] (Table 1), using the benchmarking software available as part of the SHAP package.†††https://github.com/slundberg/shap These benchmark metrics evaluate whether each attribution method finds the most important features for a given dataset and model. Expected gradients significantly outperforms the next best feature attribution method ($p=7.2\times 10^{-5}$, one-tailed Binomial test). We provide more details and additional benchmarks in Supplement Section LABEL:sec:supp_benchmark.

2.3 A pixel attribution prior improves robustness to image noise.

Prior work on interpreting image models focused on creating pixel attribution maps, which assign a value to each pixel indicating how important that pixel was for a model’s prediction [39, 47]. Attribution maps can be noisy and difficult to understand due to their tendency to highlight seemingly unimportant background pixels, indicating the model may be vulnerable to adversarial attacks [33]. Although we may prefer a model with smoother attributions, existing methods only post-process attribution maps but do not change model behavior [43, 39, 9]. Such techniques may not be faithful to the original model [13]. In this section, we describe how we applied our framework to train image models with naturally smoother attributions.

To regularize pixel-level attributions, we used the following intuition: neighboring pixels should have a similar impact on an image model’s output. To encode this intuition, we chose a total variation loss on pixel-level attributions (see Methods 3 for more detail). We applied this pixel smoothness attribution prior to the MNIST and CIFAR-10 datasets [18, 17]. On MNIST we trained a two-layer convolutional neural network; for CIFAR-10 we trained a VGG16 network from scratch (see Methods 4 for more details) [42]. In both cases we optimized hyperparameters for the baseline model without an attribution prior. To choose $\lambda$, we searched over values in $[10^{-20},10^{-1}]$ and chose the $\lambda$ that minimized the attribution prior penalty and achieved a test accuracy within $10\%$ of the baseline model. Figures 1 and 2 display expected gradients attribution maps for both the baseline and the model regularized with an attribution prior on 5 randomly selected test images on MNIST and CIFAR-10, respectively. In all examples, the attribution prior yields a model with visually smoother attributions. Remarkably, in many instances smoother attributions better highlight the target object’s structure.

Recent work has suggested that image classifiers are brittle to small domain shifts: small changes in the underlying distribution of the training and test set can significantly reduce test accuracy [30]. To simulate a domain shift, we applied Gaussian noise to images in the test set and re-evaluated the performance of the regularized and baseline models. As an adaptation of [34], we also compared the attribution prior model to regularizing the total variation of gradients with the same criteria for choosing $\lambda$. For each method, we trained 5 models with different random initializations. In Figures 1 and 2, we plot the mean and standard deviation of test accuracy on MNIST and CIFAR-10, respectively, as a function of standard deviation of added Gaussian noise. The figures show that our regularized model is more robust to noise than both the baseline and gradient-based models.

Both the robustness and more intuitive saliency maps our method provides come at the cost of reduced test set accuracy ($0.93\pm 0.002$ for the baseline vs. $0.85\pm 0.003$ for pixel attribution prior model on CIFAR-10). The trade-off between robustness and accuracy that we observe is consistent with previous work that suggests image classifiers trained solely to maximize test accuracy rely on features that are brittle and difficult to interpret [13, 49, 55]. Despite this trade-off, we find that at a stricter hyperparameter cutoff for $\lambda$ – within $1\%$ test accuracy of the baseline, rather than $10\%$ – our methods still achieve modest but significant robustness relative to the baseline. For results at different hyperparameter thresholds, more details on our training procedure, and additional experiments on MNIST, CIFAR-10 and ImageNet, see Supplement Sections LABEL:sec:supp_imagenet, LABEL:sec:ImageNetExperiments, LABEL:sec:supp_cifar and LABEL:sec:supp_mnist.

2.4 A Graph attribution prior improves anti-cancer drug response prediction.

In the image domain, our attribution prior took the form of a penalty encouraging smoothness over adjacent pixels. In other domains, there may be prior information about specific relationships between features that can be encoded as a graph (such as social networks, knowledge graphs, or protein-protein interactions). For example, prior work in bioinformatics has shown that protein-protein interaction networks contain valuable information for improving performance on biological prediction tasks [6]. Therefore, in this domain we regularized attributions to be smooth over the protein-protein feature graph analogously to the regular graph of pixels in the image.

Incorporating the $\Omega_{\textrm{graph}}$ attribution prior not only led to a model with more reasonable attributions but also improved predictive performance by letting us incorporate prior biological knowledge into the training process. We downloaded publicly available gene expression and drug response data for patients with acute myeloid leukemia (AML, a type of blood cancer) and tried to predict patients’ drug response from their gene expression [50]. For this regression task, an input sample was a patient’s gene expression profile plus a one-hot encoded vector indicating which drug was tested in that patient, while the label we tried to predict was drug response (measured by IC50, the concentration of the drug required to kill half of the patient’s tumor cells). To define the graph used by our prior, we downloaded the tissue-specific gene interaction graph for the tissue most closely related to AML in the HumanBase database [11].

A two-layer neural network trained with our graph attribution prior ($\Omega_{\textrm{graph}}$) significantly outperforms all other methods in terms of test set performance as measured by $R^{2}$ (Figure 3, see Methods for significance testing). Unsurprisingly, when we replace the biological graph from HumanBase with a randomized graph, we find that the test performance is no better than the performance of a neural network trained without any attribution prior. Extending the method proposed in [34] by applying our new graph prior as a penalty on the model’s gradients, rather than a penalty on the axiomatically correct expected gradient feature attribution, does not perform significantly better than a baseline neural network. We also observe substantially improved test performance when using the prior graph information to regularize a linear LASSO model. Finally, we note that our graph attribution prior neural network significantly outperforms graph convolutional neural networks, a recent method for utilizing graph information in deep neural networks [16].

To find out if our model’s attributions match biological intuition, we conducted Gene Set Enrichment Analysis (a modified Kolmogorov–Smirnov test). We measured whether our top genes, ranked by mean absolute feature attribution, were enriched for membership in any pathways (see Methods 5 and Supplementary Section LABEL:sec:supp_biological for more detail, including the top pathways for each model) [46]. We find that the neural network with the tissue-specific graph attribution prior captures far more biologically-relevant pathways (increased number of significant pathways after FDR correction) than a neural network without attribution priors (see Figure 3) [3]. Furthermore, the pathways our model uses more closely match biological expert knowledge, i.e., they included prognostically useful AML gene expression profiles as well as important AML-related transcription factors (see Figure 3 and Supplementary Section LABEL:sec:supp_biological) [20, 51]. These results are expected, given that neural networks trained without priors can learn a relatively sparse basis of genes that will not enrich for specific pathways (e.g. a single gene from each correlated pathway), while those trained with our graph prior will spread credit among genes. This demonstrates the graph prior’s value as an accurate and efficient way to encourage neural networks to treat functionally-related genes similarly.

2.5 A sparsity prior improves performance with limited training data.

Feature selection and sparsity are popular ways to alleviate the curse of dimensionality, facilitate interpretability, and improve generalization by building models that use a small number of input features. A straightforward way to build a sparse deep model is to apply an L1 penalty to the first layer (and possibly subsequent layers) of the network. Similarly, the Sparse Group Lasso (SGL) method penalizes all weights connected to a given feature [8, 37], while a simple existing attribution prior approach [32] penalizes the gradients of each feature in the model.

These approaches suffer from two problems. First, a feature with small gradients or first-layer weights may still strongly affect the model’s output [41]. A feature whose attribution value (e.g., integrated or expected gradient) is zero is much less likely to have any effect on predictions. Second, successfully minimizing the L1 or SGL penalty is not necessarily the best way to create a sparse model. A model that puts weight $w$ on 1 feature is penalized more than one that puts weight $\frac{w}{2p}$ on each of $p$ features. Prior work on sparse linear regression has shown that the Gini coefficient $G$ of the weights, proportional to 0.5 minus the area under the CDF of sorted values, avoids such problems and corresponds more directly to a sparse model [12, 57]. We extend this analysis to deep models by noting that the Gini coefficient can be written differentiably and used as an attribution prior.

Here, we show that the $\Omega_{\textrm{sparse}}$ attribution prior can build sparser models that perform better in settings with limited training data. We use a publicly available healthcare mortality prediction dataset of 13,000 patients [27], whose 35 features (118 after one-hot encoding) represent medical data such as a patient’s age, vital signs, and laboratory measurements. The binary outcome is survival after 10 years. Sparse models in this setting may enable accurate models to be trained with very few labeled patient samples or reduce cost by accurately risk-stratifying patients using few lab tests. We randomly sampled training and validation sets of only 100 patients each, placing all other patients in the test set, and ran each experiment 200 times with a new random sample to average out variance. We built 3-layer binary classifier neural networks regularized using L1, SGL, and sparse attribution prior penalties to predict patient survival, as well as an L1 penalty on gradients adapted for global sparsity from [34, 32]. The regularization strength was tuned from $10^{-7}$ to $10^{5}$ using the validation set for all methods (see Methods 6 and Supplement Section LABEL:sec:sparseparams).

The sparse attribution prior enables more accurate test predictions and sparser models when limited training data is available (Figure 4), with $p<10^{-4}$ and $T\geq 4.314$ by paired-samples $T$-test for all comparisons. We also plot the average cumulative importance of sorted features and find that the sparse attribution prior more effectively concentrates importance in the top few features. In particular, we observe that L1 penalizing the model’s gradients as in [32] rather than its expected gradients attributions performs poorly in terms of both sparsity and performance. A Gini penalty on gradients improves sparsity but does not outperform other baselines like SGL and L1 in ROC-AUC. Finally, we plot the average sparsity of the models (Gini coefficient) against their validation ROC-AUC across the full range of regularization strengths. The sparse attribution prior exhibits higher sparsity than other models and a smooth tradeoff between sparsity and ROC-AUC. Details and results for other penalties, including L2, dropout, and other attribution priors, are in Supplement Section LABEL:sec:supp_sparsity.

3 Discussion

The immense popularity of deep learning has driven its application in many areas with diverse, complicated domain knowledge. While it is in principle possible to hand-design network architectures to encode this knowledge, a more practical approach involves the use of attribution priors, which penalize the importance a model places on each of its input features when making predictions. Unfortunately, previous attribution priors have been limited, both theoretically and computationally. Binary penalties only specify whether features should or should not be important and fail to capture relationships among features. Approaches that only focus on a model’s input gradients change the local decision boundary but often fail to impact a model’s underlying decision-making. Attribution priors on more complicated attributions, like integrated gradients, have proven computationally difficult.

Our work advances previous work both by introducing novel, flexible attribution priors for multiple domains and by enabling the training of such priors with a newly defined feature attribution method. Our priors lead to smoother and more interpretable image models, biological predictive models that incorporate graph-based prior knowledge, and sparser healthcare models that perform better in data-scarce scenarios. Our attribution method not only enables the training of said priors, but also outperforms its predecessor – integrated gradients – in terms of reliably identifying the features models use to make predictions.

There remain many avenues for future work in this area. We chose to base our prior on an improved version of integrated gradients because it is the most prominent differentiable feature attribution method we are aware of, but a wide array of other attribution methods exist. Our framework makes it straightforward to substitute any other attribution method as long as it is differentiable, and studying the effectiveness of other attribution methods as priors would be valuable. In addition, while we develop new, more sophisticated attribution priors and show their value, there is ample room to improve on our priors and evaluate entirely new ones for other tasks. Determining the best attribution priors for particular tasks opens a further avenue of research. We believe that surveys of domain experts to establish model desiderata for particular applications will help to develop the best priors for any given situation while offering a valuable opportunity to put humans in the loop. Overall, the dual advances of sophisticated attribution priors and expected gradients enable a broader view of attribution priors: as tools to achieve domain-specific goals without sacrificing efficiency.

Methods

Methods 1 Previous attribution priors

The first instance of what we now call an attribution prior was proposed in [34], where the regularization term was modified to place a constant penalty on the gradients of undesirable features:

[TABLE]

Here, the attribution method is the gradients of the model, represented by the matrix $\frac{\partial\mathcal{L}}{\partial X}$ whose $\ell,i$th entry is the gradient of the loss at the $\ell$th sample with respect to the $i$th feature. $A$ is a binary matrix indicating which features should be penalized in which samples.

A more general interpretation of attribution priors is that any function of any feature attribution method could be used to penalize a loss function, thus encoding prior knowledge about what properties the attributions of a model should have. For some model parameters $\theta$, let $\Phi(\theta,X)$ be a feature attribution method, which is a function of $\theta$ and the data $X$. Let $\phi_{i}^{\ell}$ be the feature importance of feature $i$ in sample $\ell$. We formally define an attribution prior as a scalar-valued penalty function of the feature attributions $\Omega(\Phi(\theta,X))$, which represents a log-transformed prior probability distribution over possible attributions:

[TABLE]

where $\lambda$ is the regularization strength. Note that the attribution prior function $\Omega$ is agnostic to the attribution method $\Phi$.

Previous attribution priors [34, 19] required specifying an exact target value for the model’s attributions, but often we do not know in advance which features are important in advance. In general, there is no requirement that $\Phi(\theta,X)$ constrain attributions to particular values. Section 2 presented three newly developed attribution priors for different tasks that improve performance without requiring pre-specified attribution targets for any particular feature.

Methods 2 Expected gradients

Expected gradients is an extension of integrated gradients [47] with fewer hyperparameter choices. Like several other attribution methods, integrated gradients aims to explain the difference between a model’s current prediction and the prediction that the model would make when given a baseline input. This baseline input is meant to represent some uninformative reference input that represents not knowing the value of the input features. Although choosing such an input is necessary for several feature attribution methods [47, 41, 4], the choice is often made arbitrarily. For example, for image tasks, the image of all zeros is often chosen as a baseline, but doing so implies that black pixels will not be highlighted as important by existing feature attribution methods. In many domains, it is not clear how to choose a baseline that correctly represents a lack of information.

Our method avoids an arbitrary choice of baseline; it models not knowing the value of a feature by integrating over a dataset. For a model $f$, the integrated gradients value for feature $i$ is defined as:

[TABLE]

where $x$ is the target input and $x^{\prime}$ is baseline input. To avoid specifying $x^{\prime}$, we define the expected gradients value for feature $i$ as:

[TABLE]

where $D$ is the underlying data distribution. Since expected gradients is also a diagonal path method, it satisfies the same axioms as integrated gradients [10]. Directly integrating over the training distribution is intractable; therefore, we instead reformulate the integrals as expectations:

[TABLE]

This expectation-based formulation lends itself to a natural, sampling based approximation method: (1) draw samples of $x^{\prime}$ from the training dataset and $\alpha$ from $U(0,1)$, (2) compute the value inside the expectation for each sample and (3) average over samples.

Expected gradients also satisfies a set of important interpretability axioms: implementation invariance, sensitivity, completeness, linearity, and symmetry-preserving.

• •

Implementation invariance states that two networks with outputs that are equal over all inputs should have equivalent attributions. Any attribution method based on the gradients of a network will satisfy this axiom [47], meaning that integrated gradients, expected gradients, and gradients will all be implementation invariant.

• •

Sensitivity (sometimes called Dummy) states that when a model does not depend on a feature at all, it receives zero importance. Integrated gradients, expected gradients, and gradients all satisfy sensitivity because the gradient w.r.t. an irrelevant feature will be 0 everywhere.

• •

Completeness states that the attributions should sum to the difference between the output of a function at the input to be explained and the output of that function at a baseline. Gradients do not satisfy completeness due to saturation at the inputs; elements like ReLUs may cause gradients to be zero, making completeness impossible [47]. Integrated gradients and expected gradients both satisfy completeness due to the gradient theorem (fundamental theorem of calculus for line integrals) [47]. For expected gradients, the function being integrated is the expectation of the model’s output, so completeness means that the attributions sum to the difference between the model’s output for the input and the model’s output averaged over all possible baselines.

• •

Linearity states that for a model that is a linear combination of two submodels $f(x)=af_{1}(x)+bf_{2}(x)$, the attributions are a linear combination of the submodels’ attributions $\phi(x)=a\phi_{1}(x)+b\phi_{2}(x)$. This will hold for integrated gradients, expected gradients, and gradients because gradients are linear.

• •

Symmetry-preserving states that symmetric variables with identical values should achieve identical attributions. Integrated gradients is symmetry preserving since it is a straight line path method, and expected gradients will also be symmetry preserving, as a symmetric function of symmetric functions will itself be symmetrical [47].

The formulation of expected gradients enables an order of magnitude increase in computational efficiency relative to previous approaches for training with attribution priors. We further improve performance by reducing the need for additional data reading. Specifically, for each input in a batch of inputs, we need $k$ additional inputs to calculate expected gradients for that input batch. As long as $k$ is smaller than the batch size, we can avoid any additional data reading by re-using the same batch of input data as a reference batch, as in [56]. We accomplish this by shifting the batch of input $k$ times, such that each input in the batch uses $k$ other inputs from the batch as its reference values.

Methods 3 Specific priors

Here, we elaborate on the explicit form of the attribution priors we used in this paper. In general, minimizing the error of a model corresponds to maximizing the likelihood of the data under a generative model consisting of the learned model plus parametric noise. For example, minimizing mean squared error in a regression task corresponds to maximizing the likelihood of the data under the learned model, assuming Gaussian-distributed errors:

[TABLE]

where $\theta_{MLE}$ is the maximum-likelihood estimate of $\theta$ under the model $Y=f_{\theta}(X)+\mathcal{N}(0,\sigma)$.

An additive regularization term is equivalent to adding a multiplicative (independent) prior to yield a maximum a posteriori estimate:

[TABLE]

Here, adding an L2 penalty is equivalent to MAP for $Y=f_{\theta}(X)+\mathcal{N}(0,\sigma)$ with a $\mathcal{N}(0,\frac{1}{\lambda})$ prior. We next discuss the functional form of the attribution priors enforced by our penalties.

Methods 3.1 Pixel attribution prior

Our pixel attribution prior is based on the anisotropic total variation loss and is given as follows:

[TABLE]

where $\phi^{\ell}_{i,j}$ is the attribution for the $i,j$-th pixel in the $\ell$-th training image. Research shows [2] that this penalty is equivalent to placing 0-mean, iid, Laplace-distributed priors on the differences between adjacent pixel values, i.e., $\phi^{\ell}_{i+1,j}-\phi^{\ell}_{i,j}\sim\text{Laplace}(0,\lambda^{-1})$ and $\phi^{\ell}_{i,j+1}-\phi^{\ell}_{i,j}\sim\text{Laplace}(0,\lambda^{-1})$. [2] does not call our penalty “total variation,” but it is in fact the widely used anisotropic version of total variation and is directly implemented in Tensorflow [1, 22, 40].

Methods 3.2 Graph attribution prior

For our graph attribution prior, we used a protein-protein or gene-gene interaction network and represented these networks as a weighted, undirected graph. Formally, assume we have a weighted adjacency matrix $W\in\mathbb{R}_{+}^{p\times p}$ for an undirected graph, where the entries encode our prior belief about the pairwise similarity of the importances between two features. For a biological network, $W_{i,j}$ encodes either the probability or strength of interaction between the $i$-th and $j$-th genes (or proteins). We encouraged similarity along graph edges by penalizing the squared Euclidean distance between each pair of feature attributions in proportion to how similar we believe them to be. Using the graph Laplacian ($L_{G}=D-W$), where $D$ is the diagonal degree matrix of the weighted graph, this becomes:

[TABLE]

In this case, we choose to penalize global rather than local feature attributions. We define $\bar{\phi}_{i}$ to be the importance of feature $i$ across all samples in our dataset, where this global attribution is calculated as the average magnitude of the feature attribution across all samples: $\bar{\phi}_{i}=\frac{1}{n}\sum_{\ell=1}^{n}|\phi_{i}^{\ell}|$. Just as the image penalty is equivalent to placing a Laplace prior on adjacent pixels in a regular graph, the graph penalty $\Omega_{\textrm{graph}}$ is equivalent to placing a Gaussian prior on adjacent features in an arbitrary graph with Laplacian $L_{G}$ [2].

Methods 3.3 Sparse attribution prior

Our sparsity prior uses the Gini coefficient $G$ as a penalty, which is written:

[TABLE]

By taking exponentials of this function, we find that minimizing the sparsity regularizer is equivalent to maximizing likelihood under a prior proportional to the following:

[TABLE]

To our knowledge, this prior does not directly correspond to a named distribution. However, we observe that its maximum value occurs when one $\bar{\phi}_{i}$ is 1 and all others are [math], and that its minimum occurs when all $\bar{\phi}_{i}$ are equal. This is similar to the total variation penalty $\Omega_{\textrm{image}}$, but it is normalized and has a flipped sign to encourage differences. The corresponding attribution prior is maximized when global attributions are zero for all but one feature and minimized when attributions are uniform across features.

Methods 4 Image model experimental settings

We trained a VGG16 model from scratch modified for the CIFAR-10 dataset, containing 60,000 colored 32x32-pixel images divided into 10 categories, as in [21]. To train this network, we used stochastic gradient descent with an initial learning rate of 0.1 and an exponential decay of 0.5 applied every 20 epochs. Additionally, we used a momentum level of 0.9. For augmentation, we shifted each image horizontally and vertically by a pixel shift uniformly drawn from the range [-3, 3], and we randomly rotated each image by an angle uniformly drawn from the range [-15, 15]. We used a batch size of 128. Before training, we normalized the training dataset to have zero mean and unit variance, and standardized the test set with the mean and variance of the training set. We used $k=1$ background reference samples for our attribution prior while training. When training with attributions over images, we first normalized the per-pixel attribution maps by dividing by the standard deviation before computing the total variation; otherwise, the total variation can be made arbitrarily small without changing model predictions by scaling down the pixel attributions close to 0. See Supplement Section LABEL:sec:supp_cifar for more details.

We repeated the same experiment as above on MNIST, which contains 60,000 black-and-white 28x28-pixel images of handwritten digits. We trained a CNN with two convolutional layers and a single hidden layer. The convolutional layers each had 5x5 filters, a stride length of 1, and 32 and 64 filters total. Each convolutional layer was followed by a max pooling layer of size 2 with stride length 2. The hidden layer had 1024 units and a dropout rate of 0.5 during training [44]. Dropout was turned off when calculating the gradients with respect to the attributions. We trained with the Adam optimizer with the default parameters ($\alpha=0.001,\beta_{1}=0.9,\beta_{2}=0.999,\epsilon=10^{-8}$) [15]. We trained with an initial learning rate of 0.0001, with an exponential decay of 0.95 for every epoch, for a total of 60 epochs. For all models, we trained with a batch size of 50 images and used $k=1$ background reference sample per attribution while training. See Supplement Section LABEL:sec:supp_mnist for more details.

Methods 5 Biological experiments

Methods 5.1 Significance testing of results

To test the difference in $R^{2}$ attained by each method, we used a T-test for the means of two independent samples of scores (as implemented in SciPy) [52]. This is a two-sided test and can be applied to $R^{2}$ since $R^{2}$ is a linear transformation of mean squared error, which satisfies normality assumptions by the central limit theorem. When we compare the $R^{2}$ attained from 10 independent retrainings of the neural network to the $R^{2}$ attained from 10 independent retrainings of the attribution prior model, we find that predictive performance is significantly higher for the model with the graph attribution prior (t-statistic $=3.59$, $p=2.06\times 10^{-3}$).

To ensure that the increased performance in the attribution prior model was due to real biological information, we replaced the gene-interaction graph with a randomized graph (symmetric matrix with identical number of non-zero entries to the real graph, but entries placed in random positions). We then compared the $R^{2}$ attained from 10 independent retrainings of a neural network with no graph attribution prior to 10 independent retrainings of an neural network regularized with the random graph and found that test error was not significantly different between these two models (t-statistic $=1.25$,$p=0.23$). We also compared to graph convolutional neural networks, and found that our network with a graph attribution prior outperformed the graph convolutional neural network (t-statistic $=3.30$,$p=4.0\times 10^{-3}$).

Methods 5.2 Train/validation/test set allocation

To increase the number of samples in our dataset, we used as a feature the identity of the drug being tested, rather than one of a number of possible output tasks in a multi-task prediction. This follows from prior literature on training neural networks to predict drug response [29]. This yielded 30,816 samples (covering 218 patients and 145 anti-cancer drugs). Defining a sample as a drug and a patient, however, meant we had to choose carefully how to stratify samples into our train, validation, and test sets. While it is perfectly legitimate in general to randomly stratify samples into these sets, we wanted to specifically focus on how well our model could learn trends from gene expression data that would generalize to new patients. Therefore, we stratified samples at a patient-level rather than at the level of individual samples (e.g., no samples from any patient in the test set ever appeared in the training set). We split 20% of the total patients into a test set (6,155 samples) and then split 20% of the training data into a validation set for hyperparameter selection (4,709 samples).

Methods 5.3 Model class implementations and hyperparameters tested

LASSO. We used the scikit-learn implementation of the LASSO [48, 28]. We tested a range of $\alpha$ parameters from $10^{-9}$ to $1$, and we found that the optimal value for $\alpha$ was $10^{-2}$ by mean squared error on the validation set.

Graph LASSO. For our Graph LASSO, we used the Adam optimizer in TensorFlow [1], with a learning rate of $10^{-5}$ to optimize the following loss function:

[TABLE]

where $w\in\mathbb{R}^{d}$ is the weights vector of our linear model and $L_{G}$ is the graph Laplacian of our HumanBase network [11]. In particular, we downloaded the “Top Edges” version of the hematopoietic stem cell network, which was thresholded to only have non-zero values for pairwise interactions that had a posterior probability greater than $0.1$. We used the value of $\lambda^{\prime}$ selected as optimal in the regular LASSO model ($10^{-2}$, which corresponds to the $\alpha$ parameter in scikit-learn) and then tuned over $\nu^{\prime}$ values ranging from $10^{-3}$ to $100$. We found that a value of $10$ was optimal according to MSE on the validation set.

Neural networks. We tested a variety of hyperparameter settings and network architectures via validation set performance to choose our best neural networks, including the following feed-forward network architectures (where each element in a list denotes the size of a hidden layer): [512,256], [256,128], [256,256], and [1000,100]. We tested a range of L1 penalties on all of the weights of the network, from $10^{-7}$ to $10^{-2}$. All models attempted to optimize a least squares loss using the Adam optimizer, with learning rates again selected by hyperparameter tuning ranging from $10^{-5}$ to $10^{-3}$. Finally, we implemented an early stopping parameter of $20$ rounds to select the number of epochs of training (training was stopped after no improvement on validation error for 20 epochs, and the number of epochs was chosen based on optimal validation set error). We found that the optimal architecture (chosen by lowest validation set error) had two hidden layers of size 512 and 256, an L1 penalty on the weights of $10^{-3}$ and a learning rate of $10^{-5}$. We additionally found that 120 was the optimal number of training epochs.

Attribution prior neural networks. We next applied our attribution prior to the neural networks. First, we tuned networks to the optimal conditions described above. We then added extra epochs of fine-tuning where we ran an alternating minimization of the following objectives:

[TABLE]

[TABLE]

Following [34], we selected $\nu$ to be $100$ so that the $\Omega_{graph}$ term would initially be equal in magnitude to the least squares and L1 loss terms. We found that $5$ extra epochs of tuning were optimal by validation set error. We drew $k=10$ background samples for our attributions. To test our attribution prior using gradients as the feature attribution method (rather than expected gradients), we followed the exact same procedure, only we replaced $\bar{\phi}$ with the average magnitude of the gradients rather than the expected gradients.

Graph convolutional networks. We followed the implementation of graph convolution described in [16]. The architectures were searched as follows: in every network we first had a single graph convolutional layer (we were limited to one graph convolution layer due to memory constraints on each Nvidia GTX 1080-Ti GPU that we used), followed by two fully connected layers of sizes (512,256), (512,128), or (256,128). We tuned over a wide range of hyperparameters, including L2 penalties on the weights ranging from $10^{-5}$ to $10^{-2}$, L1 penalties on the weights ranging from $10^{-5}$ to $10^{-2}$, learning rates of $10^{-5}$ to $10^{-3}$, and dropout rates ranging from $0.2$ to $0.8$. We found the optimal hyperparameters based on validation set error were two hidden layers of size 512 and size 256, an L2 penalty on the weights of $10^{-5}$, a learning rate of $10^{-5}$, and a dropout rate of $0.6$. We again used an early stopping parameter and found that 47 epochs was the optimal number.

Methods 6 Sparsity experiments

Methods 6.1 Data description and processing

Our sparsity experiments used data from the NHANES I survey [27] and contained 35 variables (expanded to 118 features by one-hot encoding of categorical variables) gathered from 13,000 patients. The measurements included demographic information like age, sex, and BMI as well as physiological measurements like blood, urine, and vital sign measurements. The prediction task was a binary classification of whether the patient was still alive (1) or not (0) 10 years after data were gathered.

Data were mean-imputed and standardized so that each feature had 0 mean and unit variance. For each of the 200 experimental replicates, 100 train and 100 validation points were sampled uniformly at random; all other points were allocated to the test set.

Methods 6.2 Model

We trained a range of neural networks to predict survival in the NHANES data. The architecture, nonlinearities, and training rounds were all held constant at values that performed well on an unregularized network, and the type and degree of regularization were varied. All models used ReLU activations and a single output with binary cross-entropy loss; in addition, all models ran for 100 epochs with an SGD optimizer with learning rate 0.001 on the size-100 training data. The entire 100-sample training set fit in one batch. Because the training set was so small, all of its 100 samples were used for expected gradients attributions during training and evaluation, yielding $k=100$. Each model was trained on a single GPU on a desktop workstation with 4 Nvidia 1080 Ti GPUs.

**Architecture. **We considered a range of architectures, including single-hidden-layer 32-node, 128-node, and 512-node networks, two-layer [128,32] and [512,128]-node networks, and a three-layer [512,128,32]-node network; we fixed the [512,128,32] architecture for future experiments.

**Regularizers. **We tested a large array of regularizers in addition to those considered in the maintext. For details, see Supplement Section LABEL:sec:sparsemodel_supp.

Methods 6.3 Hyperparameter tuning

We selected the hyperparameters for our models based on validation performance. We searched all L1, L2, SGL and attribution prior penalties with 121 points sampled on a log scale over $[10^{-7},10^{5}]$ (Supplementary Figure LABEL:fig:supptune). Other penalties, not displayed in the maintext experiments, are discussed in Supplement Section LABEL:sec:sparseparams.

Methods 6.4 Maintext methods

**Performance and sparsity bar plots. **The performance bar graph (Figure 4, top left) was generated by plotting mean test ROC-AUC of the best model of each type (chosen by validation ROC-AUC) averaged over each of the 200 subsampled datasets, with confidence intervals given by 2 times the standard error over the 200 replicates. The sparsity bar graph (Figure 4, bottom left) was constructed using the same process, but with Gini coefficients rather than ROC-AUCs.

**Feature importance distribution plot. **The distribution of feature importances was plotted in the main text as a Lorenz curve (Figure 4, bottom right): for each model, the features were sorted by global attribution value $\bar{\phi}_{i}$, and the cumulative normalized value of the lowest $q$ features was plotted, from 0 at $q=0$ to 1 at $q=p$. A lower area under the curve indicates more features had relatively small attribution values, indicating the model was sparser. Because 200 replicates were run on small subsampled datasets, the Lorenz curve for each model was plotted using the averaged mean absolute sorted feature importances over all replicates. Thus, for a given model type, the $q=1$ point represented the mean absolute feature importance of the least important feature averaged over each replicate, $q=2$ added the mean importance for the second least important feature averaged over each replicate, and so on.

**Performance vs sparsity plot. **Validation ROC-AUC and model sparsity were calculated for each of the 121 regularization strengths and averaged over each of the 200 replicates. These were plotted on a scatterplot to show the possible range of model sparsities and ROC-AUC performances (Figure 4, top right) as well as the tradeoff between sparsity and performance.

**Statistical significance. **Statistical significance of the sparse attribution prior performance was assessed by comparing the test ROC-AUCs of the sparse attribution prior models on each of the 200 subsampled datasets to those of the other models (L1 gradients, L1 weights, SGL, and unregularized). Significance was assessed by 2-sided paired-samples $T$-test, paired by subsampled dataset. The same process was used to calculate the significance of model sparsity as measured by the Gini coefficient. Detailed tables of the resulting $p$-values and test statistics $T$ are shown in Supplement Section LABEL:sec:additionalmaintext.

Code Availability

Implementations of attribution priors for Tensorflow and PyTorch are available at https://github.com/suinleelab/attributionpriors. This repository also contains code reproducing main results from the paper.

Data Availability

The data for all experiments and figures in the paper are publicly available. The repository above contains a downloadable version of the dataset used for the sparsity experiment, as well as links to download the datasets used in the image and graph prior experiments. Data for the benchmarks was published as part of [26] and can be accessed at https://github.com/suinleelab/treeexplainer-study/tree/master/benchmark

Competing Interests

The authors declare no competing interests.

Acknowledgments

The results published here are partially based upon data generated by the Cancer Target Discovery and Development (CTD2) Network (https://ocg.cancer.gov/programs/ctd2/data-portal) established by the National Cancer Institute’s Office of Cancer Genomics.