Data-faithful Feature Attribution: Mitigating Unobservable Confounders via Instrumental Variables

Authors: qiheng  sun
haocheng  xia
jinfei  liu
Venue: NeurIPS_cc_2024_Conference
Published: 2024-05-15 19:02:49

Abstract: The state-of-the-art feature attribution methods often neglect the influence of unobservable confounders, posing a risk of misinterpretation, especially when it is crucial for the interpretation to remain faithful to the data. To counteract this, we propose a new approach, data-faithful feature attribution, which trains a confounder-free model using instrumental variables. The cluttered effects of unobservable confounders in a model trained as such are decoupled from input features, thereby aligning the output of the model with the contribution of input features to the target feature in the data generation. Furthermore, feature attribution results produced by our method are more robust when focusing on attributions from the perspective of data generation. Our experiments on both synthetic and real-world datasets demonstrate the effectiveness of our approaches.

Main paper: 
 introduction: 
 The increasing complexity and opacity of machine learning (ML) models in real-world applications boost the demand for feature attribution [8] . The feature attribution methods have been developed to help users understand why a model produces certain outputs from specific inputs. For instance, a loan applicant rejected by a bank's decision-making model might seek reasons behind the denial and what changes could potentially reverse the model's decision. Some recent studies [19, 9] have shifted the focus of feature attribution from a traditional model-centric perspective to a new perspective that is data-centric. Specifically, users may wish to assign importance values to features according to the data generation process, referring to the causal relationships through which features influence the target feature. For example, consider a medical application setting where a patient seeks to understand which personal features aggravated the illness and the cooperative impact of all features on the illness. What the patient really wants to know is how all features collaboratively contribute to the illness, rather than one diagnostic model's prediction output. These two aspects of feature attribution align with the concepts of model fidelity and data fidelity, respectively. Here, model fidelity refers to the attribution being consistent with the output of the explained model, while data fidelity pertains to the attribution being consistent with the data generation process. SHapley Additive exPlanations (SHAP) [27] and Integrated Gradients (IG) [38] are prevalent representatives of two distinct series of feature attribution methods, each uniquely satisfying critical axioms including sensitivity, implementation invariance, completeness, and symmetry [28] . These properties are essential for ensuring reasonability and fairness in feature attribution. The SHAP-based methods, grounded in game theory and particularly the Shapley value [32] , offer interpretations by evaluating all possible combinations of feature contributions in a discrete feature space. In contrast, Figure 1 : Arrows indicate direct effects. Since Ability s not directly observed and is correlated with education, the influence that should be attributed to Ability (arrow 1) is erroneously attributed to Education (arrow 2) in feature attribution. To fix this issue, Parental Education can be used as an instrumental variable for investigating the true impact of Education on Income.the IG-based methods, which are an analog of the Aumann-Shapley method [39] from cost-sharing, focus on continuous feature spaces by using path integration over gradients. It is worth noting that even when we aim for interpretations that are faithful to the data, we still rely on a model to predict the target feature values when certain features are selected or excluded in the attribution process. However, when unobservable confounders exist, both SHAP-based and IG-based methods may lead to misunderstandings if applied directly to the widely used predictive model. This is because unobservable confounders, although impacting the output, are entirely overlooked from the process of attribution. Consequently, their influence is instead attributed to other correlated features.Motivation example. As shown in Figure 1 , suppose a model predicting personal income includes Education as an input feature, and Ability serves as a confounder if the model does not incorporate it as an input feature. Because Ability has an indirect impact on Income through its influence on Education and a direct impact on Income simultaneously, the existing feature attribution methods tend to incorrectly attach the direct impact of Ability on Income to the impact of Education on Income due to their correlation. This concealed correlation may lead to incorrect attribution on the role of educational level in personal income, resulting in an overestimation of the education returns, as demonstrated in our experiments using real datasets ( § 5).In this paper, we develop a method to eliminate unobservable confounder effects in feature attribution in order to achieve a deeper understanding from the perspective of data fidelity. The instrumental variable method is widely used for causal analysis [24] . It lies in identifying features that directly affect those influenced by confounders, while not having a direct impact on the outcome themselves. By using the instrumental variable to control confounders that influence specific features, any resulting changes in the outcome variable are driven solely by how the instrumental variable alters that feature. For example, in Figure 1 , when examining the impact of Education on Income, a suitable instrumental variable could be the variable Parental Education as analyzed in appendix Section F.2. By observing the changes in Education resulting from variations in Parental Education, we can then discern the true effect of Education on Income, effectively isolating it from other confounders. Intuitively, the instrumental variable approach can help to mitigate the impact of unobservable confounders in feature attribution. However, the instrumental variable approach mainly focuses on evaluating the isolated impact of individual features influenced by unobservable confounders on the outcome variable, while feature attribution lies in considering the cooperative attribution of features. This means evaluating the combined contribution of Education and Other Features on Income.To bridge this gap, we propose using a two-stage model training with the instrumental variable that disrupts the association between confounders and other features. Specifically, the model is trained using features re-estimated through instrumental variables and collaborative variables. This ensures that the influence of confounders remains consistent despite variations in feature coalitions. Therefore, the marginal contribution of each feature, which is determined by assessing the impact on the model's output with and without the feature, is not affected by any confounders. Furthermore, the attribution value of each feature, calculated as the average of its marginal contributions across different feature coalitions, remains influenced by confounders. This alignment allows the contribution of input features to the model output to mirror their contribution in data generation to the target feature, thereby facilitating the attribution to be faithful to the data. Feature attribution involves explaining a model output by assigning attribution scores to the input instance. However, our focus is on data-faithful feature attribution, i.e., we are not trying to explain the output of a specific model but trying to explain the target feature through a model.To the best of our knowledge, we are the first to identify a crucial issue: unobservable confounders compromise feature attribution, especially when data fidelity is essential. To tackle this challenge, we propose training models free of confounders using instrumental variables, ensuring the feature attribution will remain faithful to the data. We validate the effectiveness of our proposed methods using both real and synthetic datasets, observing that our method achieves up to a 67% relative improvement over the baseline methods in terms of the error of attribution ratio metric in the real dataset. 
problem setup: 
 We aim to quantitatively assess the influence of each input feature on the target feature. This assessment can be viewed as a contribution assignment problem in the context of cooperative game theory [41] . Formally, given an explained input vector of d features x * = {x * 1 , . . . , x * d }, a baseline input x ′ , and a model f : R d → R which approximates the data generation equation for the target feature, our objective is to explain the difference in target feature, i.e., y * -y ′ , conducting data-faithful attribution for the input features. We assume x * and x ′ are of the same dimensionality d, and each entry can be either discrete or continuous. We denote by X the set of input features and Y the target feature, partitioning X into two subsets: X, which is influenced by unobserved confounders (denoted as E), and X, the set of other observable features. For clarity and convenience, we use x, y, x, x, ϵ to denote possible values within feature sets X, Y , X, X, E, respectively. For a given subset of features S, we denote the subset of the original vector of values by using S as a subscript, e.g., x S := {x i } i:i∈S .Shapley Value. Consider a set of players N = {1, . . . , d}. A coalition S is a subset of N that cooperates to complete a task. A utility function U(S) (S ⊆ N ) is the utility of a coalition S for a task. The marginal contribution of player i with respect to a coalition S is U(S ∪{i})-U(S). Shapley value is the unique metric that satisfies the properties of fair reward allocation, including balance, symmetry, additivity, and zero element [41] . It measures the expectation of marginal contribution by i in all possible coalitions. That is,SV i = 1 |N | S⊆N \{i} U(S ∪ {i}) -U(S) |N |-1 |S|.Computing the exact Shapley value requires enumerating all utilities for all player subsets. Therefore, the computational complexity of exactly calculating the Shapley value is exponential [46] .SHAP [27] utilizes the concept of Shapley values to attribute the contribution of each feature in a model. In the existing SHAP-based methods [27, 19] , the definition of utility functions for interpreting an input x can be divided into two categories [5] , condition expectation Shapley and intervention Shapley. In condition expectation Shapley, following the assumption that the features are generated according to a distribution D, the utility function is defined byU C (S) = E D [f (x)|x S = x * S ][37] based on the condition expectation of the model prediction under feature set S. In intervention Shapley, the utility function is defined byU I (S) = E D [f (x)|do(x S = x * S )] [44]where the operation do(x S = x * S ) means we intervene on the features S in variable x to be the same as the features in x * , while the features outside of S in x are influenced following the causal relationships of the features [30] . For conciseness, we omit the subscript D in the expectation term in the rest of the paper.IG is a pivotal method for model attribution [38] , particularly well-suited for deep neural networks due to its prerequisite that the model be continuously differentiable. This approach calculates the cumulative gradients along a straight-line path extending from a baseline input x ′ to the explained input x * . Mathematically, the attribution IG i assigned to a particular feature x * i for a given input x * and baseline x ′ is defined as:EQUATIONRemarkably, IG shares similarities with the Aumann-Shapley approach [39] and satisfies several essential properties, including linearity, dummy attribution, Affine Scale Invariance (ASI), proportionality, and symmetry [37] . Recent research has enhanced IG through its application to complex models and the refinement of the integrated path [1, 28] . 
misattribution with unobservable confounders: 
 In this section, we carry out a theoretical analysis to demonstrate how unobservable confounders mislead the feature attribution of both SHAP-based and IG-based methods. To show the influence of unobservable confounders in feature attribution, we first employ a simple structural equation to characterize the data-generating process, expressed asy = g( x, x) + ϵ,where g( x, x) can represent any linear or nonlinear continuous relationship involving both x and x. The equation allows us to clearly recognize the individual contributions of x and x to y, while also considering the unobserved effects encapsulated in the error term ϵ. Given that X is the set of features influenced by unobservable confounders E, it generally follows that given two data instancesx 1 and x 2 , E[ϵ| x1 ] ̸ = E[ϵ| x2 ] when x1 ̸ = x2 since X is influenced by E while E[ϵ|x 1 ] = E[ϵ|x 2 ]is valid as the feature set X is not affected by E.We discuss how unobservable confounders introduce errors in data-faithful feature attribution with a toy example of Figure 1 in condition expectation Shapley. We can assume ϵ as ability, x as the measurement of education level, x the work time in a week (i.e., the other variable), and y the weekly income. x, x, and ϵ each represents a single numerical variable. The data generation equations are defined as follows:ϵ ∼ Uniform(0, 1), x ∼ Uniform(0, 1) + ϵ, x ∼ Uniform(0, 1), y = x • x + ϵ.In this case, ability influences education, and the three features all have a direct influence on income. Consider a specific data instance x * = [ x * , x * ] = [1. 5, 1] we are curious about the contributions of the individual's education level and work hours to their income compared to the features distribution.In this example, it means assigning a value to the individual's education level x * = 1.5 and work time x * = 1 to evaluate their contribution to weekly income in comparison to the distribution of education levels and work hours of the population.First, we conduct feature attribution with a model f that is trained to fitE[y|x] = E[y| x, x] = g( x, x) + E[ϵ| x, x]which simulates the widely used supervised machine learning model training paradigm in reality. The utility derived from the model isU C (S) = E[f (x)|x S = x * S ] = E[g( x, x)|x S = x * S ] + E[ϵ|x S = x * S ]. The condition expectation Shapley value of x * is SV C x * = 1 2 {[U C ({ x * }) -U C (∅)] + [U C ({ x * , x * }) -U C ({x * })]}.Replacing the according utilities, we haveSV C x * = 1 2 {E[g( x * , x)] + E[ϵ| x * , x] -E[g( x, x)] -E[ϵ| x, x] + E[g( x * , x * )] + E[ϵ| x * , x * ] -E[g( x, x * )] -E[ϵ| x, x * ]}.Then, we conduct feature attribution for x * with the term which it really contribute to y, i.e., g( x, x). The utility derived is UC (S) = E[g( x, x)|x S = x * S ].According to the definition of condition expectation Shapley, SVC x * = 1 2 {[U C ({ x * }) -U C (∅)] + [U C ({ x * , x * }) -U C ({x * })]}.Replacing the accrrording utilities, we have SV C x * = 1 2 {E[g( x * , x)]-E[g( x, x)]+E[g( x * , x * )]-E[g( x, x * )]}.∆SV i = 1 N S⊆N \{i} |N |-1 |S| -1 {E[ϵ|x S∪{i} = x * S∪{i} ] -E[ϵ|x S = x * S ]}. Proof.Due to the limited space, please see the appendix for detailed proof. The same to the following propositions.Proposition 2. The expected error for marginal contribution of feature i in intervention Shapley with mode f trained to fitE[y|x] is E D [ϵ|do(x S∪{i} = x * S∪{i} )] -E[ϵ|do(x S = x * S )], resulting an expected deviation of attribution value by∆SV i = 1 N S⊆N \{i} |N |-1 |S| -1 {E[ϵ|do(x S∪{i} = x * S∪{i} )] -E[ϵ|do(x S = x * S )]}. Proposition 3.The expected error for attribution value of feature i using IG with model f trained tofit E[y|x] is ∆IG i = (x * i -x ′ i ) 1 α=0 ∂f (x ′ +α(x * -x ′ )) ∂xi - ∂g(x ′ +α(x * -x ′ )) ∂xi dα.Question: For each theoretical result, does the paper provide the full set of assumptions and a complete (and correct) proof? Answer: [Yes] Justification: We provide the full set of assumptions and a complete (and correct) proof. Guidelines:• The answer NA means that the paper does not include theoretical results.• All the theorems, formulas, and proofs in the paper should be numbered and crossreferenced. • All assumptions should be clearly stated or referenced in the statement of any theorems.• The proofs can either appear in the main paper or the supplemental material, but if they appear in the supplemental material, the authors are encouraged to provide a short proof sketch to provide intuition. • Inversely, any informal proof provided in the core of the paper should be complemented by formal proofs provided in appendix or supplemental material. • Theorems and Lemmas that the proof relies upon should be properly referenced. 
mitigating unobservable confounders via instrumental variables: 
 As demonstrated in Section 3, when it pertains to unobservable confounders, the prevalent feature attribution methods including SHAP and IG inevitably lead to misunderstandings that are not faithful to the data, since they rely on predictive model f trained to fit E[y|x]. This is fundamental because the trained predictive model has already associated the unobservable confounders with the input features. Therefore, it is tempting to ask: how can we decouple the confounders from their correlations with other features in the used model f ?One may think a straightforward solution is directly training a model f to fit g( x, x). Unfortunately, it is nearly impossible since we cannot remove the influence of unobservable confounders in the target feature y. To bridge the gap, we provide an alternative solution to train a model that gives the same attribution results for the input features as it is trained to fit g( x, x). Denote by ŷ = g( x, x) + E[ϵ], and f is trained to fit E[ŷ| x, x]. The influence of x and x on ŷ is identical to their impact on y, as both are encompassed within g( x, x).Example. For the toy example in Section 3, the utility calculated with f trained to fitE[ŷ| x, x] = g( x, x) + E[ϵ] is ÛC (S) = E[g( x, x)|x S = x * S ] + E[ϵ]. The condition expectation Shapley value of x * is ŜV C x * = 1 2 {[ Û({ x * }) -Û(∅)] + [ Û({ x * , x * }) -Û({x * })]}.By replacing the according utilities, we haveŜV x * = 1 2 {E[g( x * , x)] + E[ϵ] -E[g( x, x)] -E[ϵ] + E[g( x * , x * )] + E[ϵ] - E[g( x, x * )] -E[ϵ]} = 1 2 {E[g( x * , x)] -E[g( x, x)] + E[g( x * , x * )] -E[g( x, x * )]} = SV C x * .The advantage of attribution based on ŷ lies in the term E[ϵ] being constant, thereby breaking the association between ϵ and the input features.Proposition 4. From the perspective of data generation, the contributions of features in x and x to y are equivalent to their contributions to ŷ. When using the condition expectation Shapley, intervention Shapley, and Integrated Gradients (IG) methods, the attribution values of each feature in x and x are identical for both models f = g( x, x) and f = g( x, x) + E[ϵ].With Proposition 4, the problem becomes how to train the model f to fit E[ŷ| x, x] now. To achieve this, we introduce the instrumental variables.Instrumental Variable. The features that are used as instrumental variables, denoted as Ψ, can be effectively utilized in our model if they satisfy the following three key properties. 1) relevance: Ψ should correlate with X, ensuring that Ψ can serve as a reliable proxy for these features. 2) exogeneity: Ψ should be uncorrelated with the latent confounders E, ensuring that it is not influenced by these unobserved factors. 3) exclusion restriction: Ψ should influence the outcome Y solely through its effect on X. In other words, apart from its interaction with X, Ψ should not have any other direct or indirect pathways affecting Y . This ensures that the effect of Ψ on Y can be unambiguously attributed to its relationship with X. The effectiveness of IV-SHAP and IV-IG may be reduced when the three assumptions of instrumental variables are violated. However, the extent of this reduction depends on how severely the assumptions are violated.With the help of instrumental variables, we can establish the following equation by taking the expectation of y given x and ψ,E[y|x, ψ] = E[g( x, x)|x, ψ] + E[ϵ] = g( x, x) + E[ϵ]dM ( x|x, ψ),where ψ is a possible value of Ψ and dM ( x|x, ψ) is the conditional distribution of X. Given the T training data instances, the optimal parameters of model f trained to fit g( x, x) + E[ϵ] within the function space H are identified by minimizing the following objective:EQUATIONwhere L represents the loss metric we used to evaluate model performance and t is the index of specific data instance. We provide the training methods for supervised neural network models which are extensively employed in the real world in Section 4.3. Specifically, we discuss their loss functions and gradient computations when the objectives are regression and classification problems. The training steps are inspired by the two-stage training in causal effect estimation [2] and counterfactual prediction [18] . The re-estimated unconfounded values are sampled from the first-stage trained model, the sampling has little influence on the implementation and computation complexity of the second-stage model training. Therefore, the two-stage training process does not limit the method's practical usability. Due to the limited space, we provide details of model training with discrete input features and non-gradient model training in appendix Section D.Model Training for Continuous Feature Attribution. In the regression task, where the model f is a neural network trained for forecasting continuous value, it is also denoted as f θ ( x, x), where θ represents the model parameters. As our objective, we adopt a l 2 loss function. With the unknown conditional distribution of X given X and Ψ, we initially utilize a neural network model, denoted Mϕ , where ϕ is the model parameters, to approximate this distribution. The l 2 loss function for determining the optimal model parameters θ subsequently approximates as per the following equationEQUATIONwhere the integral term estimates the expected output of f θ under the distribution approximated by Mϕ . By employing the relevant calculations, we ascertain that the gradient of the loss function with respect to the t th training data point isEQUATIONIn short, our training process comprises two fundamental steps: 1) an instrumental variable method is applied to re-estimate x, mitigating the impact of unobservable confounders, and 2) this refined x is utilized to calculate the gradients for f .Model Training for Discrete Feature Attribution. In the classification task, where y is a discrete variable representing classes, we adapt the loss function to the multi-class cross-entropyEQUATIONIn this formulation, y t,r represents the true label of the t th data point in the r th category. R denotes the total number of distinct classes into which the target variable y can be classified. The probability of the model classifying a data point into the r th category is given by f θ,r ( x, x t ). The gradient calculation for the t th training data point, considering this loss function, is thenEQUATIONThis adaptation of the loss function for discrete target variables ensures that our model can handle classification tasks, effectively optimizing its performance across multiple categories.Feature Attribution Computation. The exact computation of Shapley value and integrated gradients needs huge cost while approximation methods are widely used. Towards practical applications, we further propose a Shapley value approximation method and an integrated gradients approximation method for saving computation costs in Sections E.1 and E.2, respectively. The correlation of input features may affect the data-faithfulness of IV-SHAP and IV-IG. We can combine methods which deal with the correlated input features to the two-stage model to better capture these correlations. For example, on-manifold Shapley [26] can be used to account for feature correlations, while causal Shapley [19] can be applied if the causal structure of the input features is known.Question: Does the paper fully disclose all the information needed to reproduce the main experimental results of the paper to the extent that it affects the main claims and/or conclusions of the paper (regardless of whether the code and data are provided or not)? Answer: [Yes] Justification: The experiments sections (Section 5 and F in the appendix) provided in the paper disclose all the information needed to reproduce the main experimental results.Guidelines:• The answer NA means that the paper does not include experiments.• If the paper includes experiments, a No answer to this question will not be perceived well by the reviewers: Making the paper reproducible is important, regardless of whether the code and data are provided or not. • If the contribution is a dataset and/or model, the authors should describe the steps taken to make their results reproducible or verifiable. • Depending on the contribution, reproducibility can be accomplished in various ways.For example, if the contribution is a novel architecture, describing the architecture fully might suffice, or if the contribution is a specific model and empirical evaluation, it may be necessary to either make it possible for others to replicate the model with the same dataset, or provide access to the model. In general. releasing code and data is often one good way to accomplish this, but reproducibility can also be provided via detailed instructions for how to replicate the results, access to a hosted model (e.g., in the case of a large language model), releasing of a model checkpoint, or other means that are appropriate to the research performed. • While NeurIPS does not require releasing code, the conference does require all submissions to provide some reasonable avenue for reproducibility, which may depend on the nature of the contribution. , with an open-source dataset or instructions for how to construct the dataset). (d) We recognize that reproducibility may be tricky in some cases, in which case authors are welcome to describe the particular way they provide for reproducibility.In the case of closed-source models, it may be that access to the model is limited in some way (e.g., to registered users), but it should be possible for other researchers to have some path to reproducing or verifying the results. • At submission time, to preserve anonymity, the authors should release anonymized versions (if applicable). • Providing as much information as possible in supplemental material (appended to the paper) is recommended, but including URLs to data and code is permitted. 
experiments: 
 In this section, we present our empirical evaluation in detail. We employ synthetic and real-world datasets to evaluate the faithfulness and robustness against unobservable confounders of feature attributions given by our proposed data-faithful feature attribution methods to prevalent SHAP-based and IG methods. Comparisons of feature attribution for classification problems, feature attribution on non-gradient training models, and the approximation methods of SHAP and IG are given in appendix Sections F.3, F.4, and F.5, respectively. Our code can be found in the repository at https://github.com/ZJU-DIVER/IV-SHAP.We first conducted a data simulation experiment to validate our proposed methods' effectiveness. Synthetic datasets offer an advantage in studying feature attribution, as we can obtain the ground truth of attribution values according to the data generation equation which is unobtainable in most real datasets.Dataset Generation Process. We generated two synthetic datasets, each containing four parts: an unobserved confounder ϵ, a variable x influenced by the unobserved confounder, collaborative variables x = {x i } (1 ≤ i ≤ 6), and the target feature y. Notably, dataset A and dataset B share the same x and ψ. ϵ is formulated by a uniform variable v and a parameter ρ which controls the noise level. The generation of these features adhered to specific functional relationships, as illustrated in v ∼ Uniform(0, 1),x i ∼ Uniform(0, 1) (1 ≤ i ≤ 6), ψ ∼ Uniform(0, 1), Dataset A            ϵ A = v • ρ, xA = ϵ a • ψ + ϵ A + ψ 2 /3, y A = xA + x 2 1 + x 2 + √ x 3 + x 2 4 +x5+ √ x6 2 6 + ϵ A , Dataset B                ϵ B = exp(v • ρ -1) ρ , xB = ϵ B • ψ + ϵ B + ψ 2 /3, y B = xB • exp(x 1 ) + x 2 + √ x 3 + exp(x4)+x5+ √ x6 2 6 + ϵ B .Compared Methods. We utilized two representative feature attribution algorithms, SHAP [27] and IG [38] , as the baseline methods. Specifically, we trained a neural network model on synthetic datasets to fit E[y|x] as the baseline model. Then we applied the two feature attribution methods to attribute contributions for input features. For our proposed approach, we employed the same neural network architecture but trained the model to fit E[ŷ|x] with instrumental variables. Our attribution approaches, applied to the model trained with the instrumental variable, are referred to as IV-SHAP and IV-IG, corresponding to SHAP and IG, respectively. It is worth noting that for the model, the explicitly input data features have no causal relationships among them. Therefore, Causal SHAP [19] ,Asymmetric SHAP [14] , BSHAP [37] and SHAP are equivalent in this context.Experimental Results. We randomly generated 1000 data points based on the data generation equations. We then adjusted features x and x of each data point by subtracting a certain value as a baseline input. We conducted experiments with varied subtracted values set at 0.125, 0.25, 0.375, and 0.5. For each data point, we applied IV-SHAP, IV-IG, SHAP, and IG to attribute the contributions of the features. Subsequently, we compared the attribution values of feature x against the ground truth obtained directly from the data generation equations. We observe that the errors in attribution results provided by IV-SHAP and IV-IG are significantly smaller than those of SHAP and IG. The absolute errors in the attribution values of each algorithm for every data point, as compared to the benchmark, are illustrated in Figures 2 and 3 .We conducted experiments on two real-world datasets to demonstrate the efficacy of IV-SHAP and IV-IG in practical scenarios where data generation processes are black-box. Initially, we excluded specific features to simulate unobservable confounders. Subsequently, we trained a model using the complete set of features to establish ground-truth attribution values, thereby assessing the robustness and reliability of the compared methods under realistic conditions.Real-world Datasets. The first real dataset we used is the Griliches76 dataset [17, 36] , consisting of 758 entries with 20 variables each, gathered from the U.S. labour market. This dataset is extensively used in research to explore the impact of education on income. In the study examining the relationship between the logarithm of weekly earnings (lw) and other features such as educational years (edu), years of work experience (expr), tenure at the current organization (tenure), marital status (mr), residence in the South (rns), and urban residence (smsa), there exists a significant challenge. Ability, as an unobservable confounder, not only directly influences an individual's education level but also their income. Using feature attribution methods like SHAP and IG without accounting for the confounder might incorrectly attribute the effect of ability on income to correlated educational levels.To address this issue, we incorporated the educational years of the mother (medu) as an instrumental variable to affect an individual's education. Additionally, IQ scores (iq) and knowledge in the world of work test (kww) in the dataset, serving as crucial indicators of ability, offer a unique opportunity to approximate the ground truth of the real contribution of each feature.Compared Methods. We assume a decrease in educational years for each individual as baseline inputs and execute a two-phase experiment. Initially, we omit IQ and the world of work test scores, calculating attribution values for IV-SHAP, IV-IG, SHAP, and IG. Note that the process of computing these attribution values using IV-SHAP, IV-IG, SHAP, and IG is consistent with the synthetic dataset experiments. These methods are employed to evaluate the impact of reduced education years. Subsequently, we incorporate IQ and the world of work test scores to train a new model and recalculate attribution values using SHAP and IG. Due to the absence of a real-world benchmark in the reality dataset, we adopt the attribution results from the model, which includes unobservable confounders in its training process, as our benchmark for comparison.Evaluation Metric. Denote the average attribution ratio of IV-SHAP by EAR IVSHAP =1 n n i=1 | IVSHAPi lwi |where IVSHAP i refers to the educational attribution value for the i th data point, calculated by the IV-SHAP method. IVSHAP i represents the extent to which changes in educational years influence the income in the i th data point. lw i is the income for the i th data point. EAR IVSHAP is the average of the absolute values of the ratios between the educational attribution values and income across all data points. This measure provides a comprehensive quantification of the impact of educational years on income. The average attribution ratio for reduced educational years, calculated by SHAP in the model trained with IQ and the world of work test scores, is denoted asEAR BMSHAP = 1 n n i=1 | BMSHAPi lwi| where BMSHAP i represents the benchmark attribution of education on income, incorporating IQ and the world of work test. We then compute the absolute relative error between the attributions of IV-SHAP and the benchmark using the formula | EARIVSHAP-EARBMSHAP EARBMSHAP |. For the SHAP algorithm, EAR BMSHAP still serve as a benchmark for EAR SHAP . For the attribution values calculated by IV-IG and IG, we use the average attribution ratio obtained by the IG algorithm on the model that includes IQ and kww as inputs as the benchmark.Experimental Results. The experimental results are shown in Table 2 , which demonstrates that our methods can significantly reduce the attribution error. The values in the table represent the mean of five independent runs, with the standard deviation following each mean. Empirical Analysis. The second real dataset we use is the Angrist and Krueger dataset [3] , which is an American census dataset consisting of statistical data on people born in specific years, including variables such as age (AGE), an education level (EDUC), weekly wage (LWKLYWGE), marital status (MARRIED), and race (RACE). We employed this dataset to examine the confounder effects of the ability on the attribution of education for income. In this dataset, we used the quarter of birth as an instrumental variable for years of education. The rationale behind this is the compulsory education laws in various states, which typically mandate schooling until the age of 16. Students born early in the year often start school later, leading to systematic differences in educational duration based on birth quarter. However, this dataset lacks measures of intelligence or work capability to assess the factor of ability, so we cannot conduct the experiment like the previous one. Nevertheless, the average attribution ratio of one additional year, calculated by SHAP, IV-SHAP, IG, and IV-IG, are 0.0218, 0.0206, 0.0218, and 0.0207, respectively. This aligns with the observation that people may overestimate educational returns because of neglecting the confounder ability in [6] . 
limitations: 
 Despite the strengths of our approach, there are several limitations to consider which are shown as follows:• Dependence on the Availability of Instrumental Variables: Our approach assumes the presence of suitable instrumental variables for features affected by unobserved confounders. However, in practical scenarios, finding appropriate instrumental variables can be challenging sometimes. For further information on identifying instrumental variables, refer to works such as [4] , [10] , and [23] .• Linearity Assumption in Theoretical Derivations: Our theoretical derivations are based on the assumption that the influence of unobserved confounders on the target features is linear. This assumption does not hold in all real-world situations. Nevertheless, in our experiments with real datasets in Section 5.2, our attribution method showed significant improvements over existing methods, even when the influence of unobservable confounders on features was non-linear.These limitations highlight areas for future research, particularly in developing methods that do not rely on the availability of instrumental variables and that can give theoretical analysis of nonlinear effects of unobserved confounders. Addressing these aspects can enhance the practicality and applicability of our methods.Question: Does the paper specify all the training and test details (e.g., data splits, hyperparameters, how they were chosen, type of optimizer, etc.) necessary to understand the results?Answer: [Yes]Justification: The full details are provided with the code.Guidelines:• The answer NA means that the paper does not include experiments.• The experimental setting should be presented in the core of the paper to a level of detail that is necessary to appreciate the results and make sense of them. • The full details can be provided either with the code, in appendix, or as supplemental material. 
conclusion: 
 In this paper, we focus on addressing the effects of unobservable confounders in feature attribution, emphasizing feature attribution being faithful to data. The proposed method improves the understanding of the causal factors driving an outcome variable, going beyond standard attribution scores that simply describe predictive models. Our approach of training confounder-free models using instrumental variables effectively isolates the impact of confounders, enhancing the robustness of data-faithful feature attribution results. Our validations using real and synthetic datasets confirm the effectiveness of the proposed methods. For future work, we intend to develop methods that do not rely on the availability of instrumental variables and that can provide a theoretical analysis of the non-linear effects of unobserved confounders. For the broader impacts of the paper, please see Section A in the appendix due to the limited space.Question: Does the paper report error bars suitably and correctly defined or other appropriate information about the statistical significance of the experiments?Answer: [Yes]Justification: For the results on synthetic datasets, we used bar charts, which effectively demonstrate the stability of the experimental results. For the experimental results presented in tabular format, we included the standard deviation to show statistical significance.Guidelines:• The answer NA means that the paper does not include experiments.• The authors should answer "Yes" if the results are accompanied by error bars, confidence intervals, or statistical significance tests, at least for the experiments that support the main claims of the paper. • The factors of variability that the error bars are capturing should be clearly stated (for example, train/test split, initialization, random drawing of some parameter, or overall run with given experimental conditions). • The method for calculating the error bars should be explained (closed form formula, call to a library function, bootstrap, etc.) • The assumptions made should be given (e.g., Normally distributed errors). • It should be clear whether the error bar is the standard deviation or the standard error of the mean. • It is OK to report 1-sigma error bars, but one should state it. The authors should preferably report a 2-sigma error bar than state that they have a 96% CI, if the hypothesis of Normality of errors is not verified. • For asymmetric distributions, the authors should be careful not to show in tables or figures symmetric error bars that would yield results that are out of range (e.g. negative error rates). • If error bars are reported in tables or plots, The authors should explain in the text how they were calculated and reference the corresponding figures or tables in the text. 
experiments compute resources: 
 Question: For each experiment, does the paper provide sufficient information on the computer resources (type of compute workers, memory, time of execution) needed to reproduce the experiments?Answer: [Yes] Justification: Please see Section F.1 in the appendix for the information of computer resources.Guidelines:• The answer NA means that the paper does not include experiments.• The paper should indicate the type of compute workers CPU or GPU, internal cluster, or cloud provider, including relevant memory and storage. • The paper should provide the amount of compute required for each of the individual experimental runs as well as estimate the total compute. • The paper should disclose whether the full research project required more compute than the experiments reported in the paper (e.g., preliminary or failed experiments that didn't make it into the paper). Guidelines:• The answer NA means that there is no societal impact of the work performed.• If the authors answer NA or No, they should explain why their work has no societal impact or why the paper does not address societal impact. • Examples of negative societal impacts include potential malicious or unintended uses (e.g., disinformation, generating fake profiles, surveillance), fairness considerations (e.g., deployment of technologies that could make decisions that unfairly impact specific groups), privacy considerations, and security considerations. • The conference expects that many papers will be foundational research and not tied to particular applications, let alone deployments. However, if there is a direct path to any negative applications, the authors should point it out. For example, it is legitimate to point out that an improvement in the quality of generative models could be used to generate deepfakes for disinformation. On the other hand, it is not needed to point out that a generic algorithm for optimizing neural networks could enable people to train models that generate Deepfakes faster. • The authors should consider possible harms that could arise when the technology is being used as intended and functioning correctly, harms that could arise when the technology is being used as intended but gives incorrect results, and harms following from (intentional or unintentional) misuse of the technology. • If there are negative societal impacts, the authors could also discuss possible mitigation strategies (e.g., gated release of models, providing defenses in addition to attacks, mechanisms for monitoring misuse, mechanisms to monitor how a system learns from feedback over time, improving the efficiency and accessibility of ML). 
safeguards: 
 Question: Does the paper describe safeguards that have been put in place for responsible release of data or models that have a high risk for misuse (e.g., pretrained language models, image generators, or scraped datasets)?14. Crowdsourcing and Research with Human Subjects Question: For crowdsourcing experiments and research with human subjects, does the paper include the full text of instructions given to participants and screenshots, if applicable, as well as details about compensation (if any)? Answer:[NA] Justification: The paper does not involve crowdsourcing nor research with human subjects. We used two real-world datasets that contain information about education and income collected from the market, which are anonymized and do not involve direct human subjects research. Therefore, the requirements for providing instructions, screenshots, and compensation details do not apply. Guidelines:• The answer NA means that the paper does not involve crowdsourcing nor research with human subjects. • Including this information in the supplemental material is fine, but if the main contribution of the paper involves human subjects, then as much detail as possible should be included in the main paper. • According to the NeurIPS Code of Ethics, workers involved in data collection, curation, or other labor should be paid at least the minimum wage in the country of the data collector. 
institutional review board (irb) approvals or equivalent for research with human subjects: 
 Question: Does the paper describe potential risks incurred by study participants, whether such risks were disclosed to the subjects, and whether Institutional Review Board (IRB) approvals (or an equivalent approval/review based on the requirements of your country or institution) were obtained? Answer: [NA] Justification: The paper uses two real-world datasets that contain information about education and income collected from the market. These datasets are anonymized and do not involve direct human subjects research. Therefore, IRB approval is not required. Guidelines:• The answer NA means that the paper does not involve crowdsourcing nor research with human subjects. • Depending on the country in which research is conducted, IRB approval (or equivalent) may be required for any human subjects research. If you obtained IRB approval, you should clearly state this in the paper. • We recognize that the procedures for this may vary significantly between institutions and locations, and we expect authors to adhere to the NeurIPS Code of Ethics and the guidelines for their institution. • For initial submissions, do not include any information that would break anonymity (if applicable), such as the institution conducting the review. 
 
 Appendix: In the appendix of our paper, we provide comprehensive additional content. We discuss the lbroader impacts of the paper in Section A, respectively. Section B reviews the related works. In Section D, we delve into the computation of gradients when dealing with discrete features and discuss training methods for non-neural network models. Following this, Section E presents the algorithm for computing Shapley values optimized using confidence intervals, along with an error analysis for unbiased sampling in integrated gradients. Subsequently, Section F offers supplementary material related to our experimental procedures. This includes a detailed analysis of the characteristics of our experimental dataset, justifying our experimental design. Additional results from classification experiments and non-neural network models are also provided.While we believe our paper has many positive social impacts, we think it can particularly affect:• Fairness and Equity in Automated Systems: Reduces biases caused by unobservable confounders in feature attribution, promoting fairness in systems like credit scoring and hiring. This helps to build trust in these systems and supports fair decision-making in various areas.• Improved Decision-Making in Healthcare: Our method enables more accurate identification of factors affecting patient outcomes, leading to better diagnosis, treatment plans, and personalized medicine. Healthcare professionals can make better decisions, which improves patient care and outcomes.We do not think our paper has any negative social impacts.In this section, we first introduce seminal works that have a significant influence on the feature attribution domain. This is followed by an exploration of works integrating causal knowledge into feature attribution. Additionally, we introduce the widespread presence of confounders in machine learning. Finally, we discuss advancements in computational optimization for SHAP-based methods.For a more detailed survey of feature attribution, please see [12] .LIME (Local Interpretable Model-agnostic Explanations) facilitates the understanding of individual predictions of complex models by creating explanatory models [31] . It reveals the impact of features on predictions by perturbing the input and observing the resultant changes in output. DeepLIFT (Deep Learning Important FeaTures) offers a method for assessing feature importance in deep neural networks by comparing the activation of each feature against a reference activation, proving particularly effective in interpreting deep learning models [33] . SmoothGrad enhances visual interpretations of gradient-based methods by applying multiple small random perturbations to the input data and averaging the gradients of these perturbations [35] . Meanwhile, researchers have increasingly recognized that for interpretability methods to be effective and credible, they need to satisfy axiomatic properties. SHAP (SHapley Additive exPlanations) employs Shapley values from cooperative game theory to measure feature contributions, offering a model-agnostic approach with broad applicability [27] . It adheres to desirable allocation properties, ensuring both consistency and equity in attributing feature influence on predictions. Meanwhile, Integrated Gradients (IG) calculates feature importance through the integration of gradients along a straight path from a baseline to the input, making it ideal for scenarios with continuous features and differentiable models [38] .In the evolving field of feature attribution, the significance of causal relationships for data-faithful interpretations is increasingly recognized [21] . Asymmetric Shapley values (ASVs) are developed to infuse causal understanding into model explanations [14] . They achieve this by modifying the symmetry axiom in the Shapley value framework, allowing for the inclusion of causal relationships. Notably, ASVs can provide insights even without a complete causal graph. Causal Shapley values stand out in their capacity to distinguish between direct and indirect feature impacts on model predictions, offering a profound understanding of data generation [19] . Shapley Flow distinguishes itself by evaluating the entire causal graph, attributing influence across its edges rather than focusing solely on nodes [40] . Recursive Shapley Value (RSV) presents a specialized approach for graphical models, quantifying the propagation of changes from source nodes throughout the graph [34] . However, despite the advancements made by these methods in considering the causal relationships between model input features, these methods overlook the impact of unobservable confounders on feature attribution.Researchers have recently begun exploring methods to identify and adjust for confounders in algorithmic models to enhance decision-making quality [22, 42] . Gao et al. [15] identify that pre-trained graph neural networks perform better on pruned graphs than on full graphs due to confounders and introduce Robust Causal Graph Representation Learning (RCGRL) to effectively address this issue by eliminating confounders. Zhang et al. [45] introduce a method for causal imitation learning in the presence of unobservable confounders, featuring a graphical criterion to evaluate its feasibility despite partially observed decision variables behind expert actions. Deep Sequential Weighting (DSW) is proposed for estimating individual treatment effects in healthcare, accounting for time-varying hidden confounders using deep learning [25] . In confounded sequential decision-making, Xu et al. [43] study introduces an instrumental variable (IV) method for off-policy evaluation (OPE) to estimate policy returns accurately in infinite horizon settings.TreeSHAP [26] , for tree-based models, enhances SHAP value calculation efficiency by utilizing tree structures to skip redundant feature combination evaluations. Dynamic Shapley, as discussed in the paper by [47] , focuses on dealing with scenarios where the players may change. Kernel SHAP [27] , suitable for various models, approximates SHAP values by sampling in the feature space and assessing the impact of different feature combinations. Among the recent advancements in optimizing SHAP computation are TMC (Truncated Monte Carlo) [16] and FastSHAP [20] , each offering unique approaches to enhance efficiency. TMC employs a truncation technique for rapid, biased sampling approximations. FastSHAP is biased too, employing a pre-trained auxiliary model, speeds up SHAP value prediction. Moreover, unlike methods that approximate Shapley values through sampling, its ability to accurately estimate Shapley values does not improve with more samples, as the precision of the auxiliary model is predetermined upon training. Recently, researchers have proposed one Shapley value approximation method based on the complementary contribution which can be adapted to the general class of feature attribution scenarios [46] .Proof. From the perspective of the data generation process for y, the marginal contribution of feature i is exclusively linked to the function g. Thus, the marginal contribution of feature i in condition expectation Shapley for the target feature generation is UC (S ∪ {i}) -U C (S) = E[g( x, x)|x S∪{i} = x S∪{i} ] -E[g( x, x)|x S = x S ]. For the model trained to fit E[y| x, x], we have E[y| x, x] = E[g( x, x)| x, x] + E[ϵ| x, x] = g( x, x) + E[ϵ| x].Therefore, given explained input x, the marginal contribution of a particular feature i with S in condition expectation Shapley derived with model f can be represented as followsU C (S ∪{i})-U C (S) = E[f (x)|x S∪{i} = x S∪{i} ]-E[f (x)|x S = x S ] = E[g( x, x)|x S∪{i} = x S∪{i} ] + E[e|x S∪{i} = x S∪{i} ] -E[g( x, x)|x S = x S ] -E[e|x S = x S ].Thus, the marginal contribution calculated by the model f which is trained to fit E[y| x, x] includes an error termE D [ϵ|x S∪{i} = x S∪{i} ] -E[ϵ|x S∪{i} = x S∪{i} ],arises from the model's reliance on the correlations within features X and E. By averaging the errors in the marginal contributions of feature i with all possible cooperate coalitions, we can get the expected deviation of attribution value by∆SV i = 1 N S⊆N \{i} |N |-1 |S| -1 {E D [ϵ|x S∪{i} = x * S∪{i} ] -E[ϵ|x S = x * S ]}.Proof. From the data generation perspective, the marginal contribution of feature i in intervention Shapley for the target feature generation should be UI (S ∪ {i}) -U I (S) = E[g( x, x)|do(x S∪{i} = x S∪{i} )] -E[g( x, x)|do(x S = x S )].However, for the intervention Shapley, the marginal contribution derived with a model trained to fitE[y| x, x] is U I (S ∪ {i}) -U I (S) = E[f (x)|do(x S∪{i} = x S∪{i} )] -E[f (x)|do(x S = x S )] = E[g( x, x)|do(x S∪{i} = x S∪{i} )] + E[ϵ|do(x S∪{i} = x S∪{i} )] -E[g( x, x)|do(x S = x S )] -E D [ϵ|do(x S = x S )]. Thus, we have the expected error for marginal contribution of feature i in intervention Shapley with mode f trained to fitE[y|x] is E D [ϵ|do(x S∪{i} = x * S∪{i} )] -E D [ϵ|do(x S = x * S )].By averaging the errors in all the marginal contributions of feature i with all possible cooperative coalitions, we can get the expected deviation of attribution value by∆SV i = 1 N S⊆N \{i} |N |-1 |S| -1 {E D [ϵ|do(x S∪{i} = x * S∪{i} )] -E D [ϵ|do(x S = x * S )]}.Proof. For IG, where ϵ is the unobservable confounder correlated with x, the derivative ∂f ∂x is likely to incorporate the effect of ϵ on x, as it is trained to fit g(x) + ϵ. Consequently, the following inequality typically holds ∂f ∂x ̸ = ∂g ∂x . Therefore, the attribution IG i (x, x ′ , f ) derived from the predictive model f generally differs from the attribution IG i (x, x ′ , g) that should be obtained based on the actual data generation process. When we accumulate the difference of ∂f ∂x and ∂g ∂x in the path, we can get the error for attribution value of feature i using IG with model f trained to fit E[y|x] is∆IG i = (x * i -x ′ i ) 1 α=0 ∂f (x ′ +α(x * -x ′ )) ∂xi - ∂g(x ′ +α(x * -x ′ )) ∂xi dα.Proof. The marginal contribution of a particular feature i with S in condition expectation Shapley derived with model f = g( x, x) + E[ϵ] can be represented as followsU C (S ∪ {i}) -U C (S) = E[f (x)|x S∪{i} = x S∪{i} ] -E[f (x)|x S = x S ] = E[g( x, x)|x S∪{i} = x S∪{i} ] + E[e] - E[g( x, x)|x S = x S ] -E[e] = E[g( x, x)|x S∪{i} = x S∪{i} ] -E[g( x, x)|x S = x S ]. Thus, the attribution are identical for models f = g( x, x) and f = g( x, x) + E[ϵ] in condition expectation Shapley.The marginal contribution derived with a model trained to fit f= g( x, x) + E[ϵ] is U I (S ∪ {i}) - U I (S) = E[f (x)|do(x S∪{i} = x S∪{i} )] -E[f (x)|do(x S = x S )] = E[g( x, x)|do(x S∪{i} = x S∪{i} )] + E[ϵ] -E[g( x, x)|do(x S = x S )] -E[ϵ] = E[g( x, x)|do(x S∪{i} = x S∪{i} )] - E[g( x, x)|do(x S = x S )].Thus, the attribution are identical for models f = g( x, x) andf = g( x, x) + E[ϵ] in intervention Shapley.The derivative ∂f ∂x = ∂g ∂x holds when model f = g(x)+E[ϵ] as E[ϵ] is a constant. Thus, the attribution are identical for models f = g( x, x) and f = g( x, x) + E[ϵ] in Integrated Gradients(IG).We extend our discussion to scenarios where the targets predicted by Mϕ ( x | x t , ψ t ) are discrete. In cases where the prediction features P of Mϕ are discrete, the fundamental approach to optimizing the loss function L(T ; θ) remains similar. The primary modification involves substituting the integral over the probability distribution of P with a summation across discrete points. Assuming P has K categories, and denoting the probability of the t th data point being classified into the k th category by Mϕ xk | x t , z t , the loss function when y is continuous is reformulated as:EQUATION)For the t th training data point, the gradient of this loss function is:∇ θ L t = -2 y t - K k=1 Mϕ xk | x t , ψ t f θ xk , x t • K k=1 Mϕ xk | x t , ψ t f ′ θ xk x t . (8)Furthermore, the gradients of a mini-batch comprising m training data tuples are computed as:EQUATIONIn situations where y is a discrete variable, representing categories or classes, the multi-class crossentropy can be formulated as:EQUATIONIn this formulation, y t,r represents the true label of the t th data point in the r th category. R denotes the total number of distinct categories or classes into which the target variable y can be classified. The model's prediction for this category is given by f θ,r ( x, x t ). The gradient calculation for the t th training data point, considering this loss function, is then:EQUATIONFurthermore, the gradients of a mini-batch comprising m training data tuples are computed as:EQUATION)This adaptation of the loss function for discrete target variables ensures that our model can handle classification tasks, effectively optimizing its performance across multiple categories.When training models to fit ŷ in scenarios where gradient-based optimization is not feasible, we introduce an alternative approach that effectively addresses the influence of confounding factors. The essence of this approach lies in the generation of synthetic data, which is derived from the predicted distribution of p. By sampling each original data point B times, we create B synthetic data points for every original point. This process results in a synthetic dataset that embodies the controlled effects of the confounders. The creation of this dataset is a vital step towards ensuring that the subsequent model training is less influenced by confounding variables.A key advantage of this method is its independence from any specific model type. The generated synthetic dataset can be utilized to train a variety of machine learning models, not limited to those that rely on gradient-based optimization. This model-agnostic nature significantly widens the applicability of our approach, making it suitable for various scenarios and models. Through this method, we ensure that the training of models occurs in an environment where the impact of confounders is mitigated, thereby enhancing the reliability of the feature attribution.SHAP-based methods and IG-based methods can be applied to the proposed confounder-free models. However, the computational complexity poses a significant barrier to real-world applications. The exact computation of the Shapley value is proved to be an #P-hard problem [11] , and the exact computation of integrated gradients requires the antiderivative of the gradient, which is infeasible due to their complexity, necessitating the use of approximation methods. The approximation cost of the Shapley value is higher than the straightforward sampling in the path of integrated gradients due to extensive feature subset evaluations. To further enhance the applicability, we develop optimizations for the approximation of SHAP-based methods. Research has shown that Shapley values can be represented not only based on marginal contributions but also complementary contributions, which allows for reusing samples in estimations, offering an advantage [46] . We propose an enhanced approach for SHAP-based methods, optimizing complementary contribution-based sampling using confidence intervals. This optimization is designed to minimize estimation errors in all utility functions within SHAP-based methods.Recent work [46] suggests that the Shapley value formula can be equivalently transformed into a form expressed based on complementary contributions. Here, the complementary contribution refers to the difference in utility between complementary subsets. The Shapley expression is given byEQUATION)The formulas based on complementary contributions offer advantages in terms of sample reusability during approximate computation. Building on this foundation, we propose a dynamic sampling adjustment based on the confidence intervals of Shapley value estimates in stratified sampling.Denote by S i,j N = {S ∪ {z i }|S ⊆ N \ {z i }, |S| = j -1} (1 ≤ j ≤ n) the set of (z i , j)-coalitions, and by SV i,j the expected complementary contributions of (z i , j)-coalitions. That is,EQUATION)Complementary contributions CC(S) = U(S ∪ {z i }) -U(N \ (S ∪ {z i })) are naturally stratified into n strata S i,1 N , . . . , S i,n N according to the coalition size. We start by deriving the confidence interval of the estimator of SV i,j using t-test. Lemma 5. According to the Central Limit Theorem, the sample mean approximates a normal distribution when the sample size is sufficiently large. Assuming a confidence level of α, the confidence interval for SV i,j based on the t-test is SV i,j ± A α Si,j √ mi,j , which can also be represented asEQUATION)where A α is the t-score corresponding to α, m i,j is the sample size of S i,j N and S i,j is the sampling variance of SV i,j .Denote by S j N = {S|S ⊆ N , |S| = j} the set of j-coalitions (1 ≤ j ≤ n). After drawing a coalition S from S j N , we can estimate the complementary contribution CC N (S), which can be used in SV i,j for z i in S and SV i,n-j for z i in N \ S. In light of the fact that each sample can influence multiple strata, how should we allocate the number of samples to optimize the precision of the estimated values? We denote I j N as the sum of the confidence intervals for strata which can be influenced by a random sample of S j N , that isEQUATION)For each sampling iteration, we select the stratum that maximizes the sum of the corresponding confidence intervals.Next, we will outline the algorithmic procedure, which is divided into two primary stages. Due to page limitations, the pseudo-code of the algorithm is presented in Algorithm 1 in the appendix. In the first stage, we sample at least m init samples for SV i,j . We then compute unbiased estimations of σ 2 i,j using Bessel's correction based on samples collected in the first stage. In the second stage, the stratum for each individual sampling is determined based on the sum of confidence intervals I j N (1 ≤ j ≤ n/2). Furthermore, this sum is updated upon the completion of each sampling. Specifically, let CC N (S 1 ∪ {z i }), . . . , CC N (S mi,j ∪ {z i }) be m i,j samples for computing SV i,j , then σ 2 i,j =1 mi,j -1 mi,j k=1 (CC N (S k ∪ {z i }) -1 mi,j mi,j k=1 CC N (S k ∪ {z i })). Let m first be the number of samples used in the first stage, and the number of remaining samples is m -m first . We calculate I j N (1 ≤ j ≤ n/2) according to Equation ( 17) using the unbiased sample variance σ 2 i,j . We randomly sample from the stratum with the largest sum of confidence intervals. Following this sampling, the sum of confidence intervals will be updated. We will continue to repeat this process until all samples have been utilized. The final estimation of Shapley value is the average of all complementary contribution means in each stratum.SHAP experiences an exponential increase in computational complexity with the addition of more features. However, for IG, the increase in features does not significantly escalate the complexity of the integral path. Therefore, in terms of efficiency in approximate computations, IG generally outperforms SHAP. This makes IG a preferable choice in situations where computational complexity for interpretability is a critical concern. However, it is important to note that SHAP, as a modelagnostic method, is applicable for interpreting models that are not based on gradient optimization.Algorithm of SHAP Computation Based on the Confidence Interval. For the algorithm process we propose, which utilizes confidence intervals to optimize the calculation of Shapley values, refer to Algorithm 1. It is important to note that our algorithm is utility function-agnostic [29] , meaning it can be applied across various SHAP-based algorithm variants. This is achieved by simply substituting the utility function defined by each method into our calculation. Furthermore, we provide an unbiased proof of our method in Theorem 6. Theorem 6. Given a set of players N = {z 1 , . . . , z n }, Algorithm 1 gives an unbiased estimation of Shapley value for every player, that is,E[SV i ] = SV i (1 ≤ i ≤ n).Proof. Denote by CC N (S 1 ), . . . , CC N (S mi,j ) a sample of S i,j N (1 ≤ i, j ≤ n) drawn by Algorithm 1. The expectation of the sample SV i,j = 1 mi,j mi,j k=1 CC N (S k ). We can compute the expectation of SV i,j withEQUATIONAccording to Equation 15,E[CC N (S k )] = SV i,j . Thus, E[SV i,j ] = SV i,j that means SV i,j is an unbiased estimation of SV i,j .Then, we can compute the expectation of SV i produced by Algorithm 1. We haveEQUATIONThat is, SV i is an unbiased estimation of SV i .In existing literature related to Integrated Gradients (IG), interpolation methods are commonly used for approximation [13] , which already exhibit high efficiency compared to SHAP-like approximation methods. In contrast, our paper introduces the use of Monte Carlo methods for the integration of gradients. It's important to clarify that the aim of proposing an unbiased estimate for IG is not Algorithm 1 Shapley value computation based on the confidence interval.Input: players N = {z 1 , . . . , z n }, m init > 1 , and m > 0 Output: approximate Shapley value SV i for each player z i (1≤ i ≤ n) SV i , SV i,j , m i,j ← 0 (1 ≤ i ≤ n); c ← -1; while c ̸ = n j=1 m 1,j do c = n j=1 m 1,j ; for i=1 to n, j = 1 to n do if m i,j < m init then let S be a sample drawn from S j N ; u ← U(S) -U(N \ S); for z i ∈ S do SV i,|S| + = u; m i,|S| + = 1; end for for z i ∈ N \ S do SV i,|N \S| -= u; m i,|N \S| + = 1; end for end if end for end while compute S 2 i,j (1 ≤ i, j ≤ n); m f irst ← n j=1 m 1,j ; for k = 0 to m -m f irst do for j=1 to n do I j N = 2 * ( N i=1 A α Si,j √ mi,j + N i=1 A α Si,n-j√ mi,n-j ); end for let S be a sample drawn from S j N where j corresponding to the stratum with the maximumI j N ; u ← U(S) -U(N \ S); for z i ∈ S do SV i,|S| + = u; m i,|S| + = 1; end for for z i ∈ N \ S do SV i,|N \S| -= u; m i,|N \S| + = 1; end for update S 2 i,j (1 ≤ i, j ≤ n); end for for i=1 to n do SV i = 1 n n j=1 SV i,j /m i,j ; end for return SV 1 , . . . , SV n .to optimize sampling efficiency but rather to provide an error analysis for the sampling process. Our proposed approach not only provides an unbiased estimate of the integrated gradients but also includes this crucial error analysis.Let u be a uniformly distributed random variable over the interval [0, 1], whereh(u) = (x i -x ′ i ) ∂f (x ′ + u (x -x ′ )) ∂x i .Then, h(u) is an unbiased estimator ofEQUATIONNext, we can obtain a more accurate unbiased estimate through the following steps: First, we generatem i random samples u 1 , • • • , u mi of u.Then, we sample h(u) based on these random numbers to get an independent and identically distributed sample h(u 1 ), • • • , h(u mi ). Finally, we use the observed values ofIG i (x, x ′ , f ) = 1 mi mi i=1 h(u i ) as the estimate for IG i (x, x ′ , f ). The proof that 1 mi mi i=1 h(u i ) is an unbiased estimator of IG i (x, x ′ , f ) is straightforward, due to the fact that E[ 1 mi mi i=1 h(u i )] = 1 mi mi i=1 E[h(u i )] = IG i (x, x ′ , f ). Lemma 7. The probability that IG i (x, x ′ , f )(1 ≤ i ≤ n) deviates from IG i (x, x ′ , f) be equal to or greater than any fixed ϵ ≥ 0 given the sample size m i is bounded byEQUATIONwherer i,j = max u∈[0,1] h(u) -min u∈[0,1] h(u).Proof. According to Hoeffding's inequality, we haveEQUATIONEQUATION=P(| mi i=1 h(u i ) -E[ mi i=1 h(u i )]| ≥ m i ϵ|m i ) ≤ 2 exp(- 2m i ϵ 2 r 2 i,j).F Supplement to ExperimentsWe conduct experiments on a machine with 2 Montage(R) Jintide(R) C6226R @ 2.90GHz and 256GB memory. Our experiments do not require high-end hardware, and our algorithm is not time-consuming. For feature attribution experiments on synthetic datasets, each attribution algorithm takes about 10 seconds to attribute one data point. Thus, an attribution algorithm takes several hours to attribute an entire dataset. The time consumption on our real dataset is similar. Also, since our algorithm uses very little memory, it is easy to run multiple processes and algorithms in parallel on a machine, so the time cost for reproduction is friendly.Statistical Characteristics of the Synthetic Datasets. We present the mean and variance of each feature in our synthetic datasets to provide crucial insights into their characteristics, as shown in Table 3 . The mean offers an understanding of the average behaviour of features, while the variance indicates their variability. This information is vital for assessing the data's overall distribution and quality, and it plays a key role in interpreting the results of our proposed methods and baseline algorithms. .0 e a 0.10 ± 0.05 0.20 ± 0.11 0.30 ± 0.17 0.40 ± 0.23 0.50 ± 0.28 t a 0.21 ± 0.12 0.27 ± 0.14 0.32 ± 0.17 0.37 ± 0.18 0.42 ± 0.21 y a 0.68 ± 0.18 0.87 ± 0.24 1.0 ± 0.31 1.15 ± 0.39 1.30 ± 0.47 e b 2.03 ± 0.12 1.12 ± 0.13 0.83 ± 0.15 0.69 ± 0.16 0.61 ± 0.18 t b 1.11 ± 0.21 0.72 ± 0.19 0.60 ± 0.18 0.53 ± 0.18 0.50 ± 0.18 y b 2.83 ± 0.25 1.65 ± 0.23 1.27 ± 0.24 1.09 ± 0.25 0.99 ± 0.27 These statistics help validate the synthetic data's consistency and reliability, which is essential for the credibility of our experimental findings.Statistical Characteristics of the Real Datasets. In the analysis of the Griliches76 dataset [36] , we observe various degrees of correlation between the mother's years of education (denoted as med) and several key variables. Firstly, the correlation coefficient between the mother's education and marital status (variables mrt and mrt80) is close to 0, indicating almost no correlation between the mother's level of education and her marital status. The correlation coefficients with urban residence status (variables smsa and smsa80) are 0.098 and 0.031, respectively, suggesting a slight positive correlation. This implies that there is a weak but positive association between the mother's educational attainment and living in an urban area. A moderate positive correlation is observed with IQ, as indicated by a correlation coefficient of 0.226 with the mother's education. This suggests that higher maternal education is somewhat associated with higher IQ scores. Similarly, the correlation between the mother's education and scores in the world of work test is 0.195, which also reflects a moderate positive correlation. This indicates that higher maternal education levels might be linked to better performance in job-related knowledge. Most notably, the correlation coefficients with personal education years (variables s and s80) are 0.340 and 0.341, respectively, indicating a relatively strong positive correlation. This suggests that the mother's level of education is considerably associated with the individual's own educational attainment. While there is a certain degree of correlation between maternal education and both IQ and job knowledge test scores, factors like regression to the mean in intelligence suggest that the correlation between a mother's education and her child's IQ is weaker than the correlation between a mother's education and the child's own educational attainment. Therefore, we posit that selecting maternal education as an instrumental variable, although not perfectly ideal, still holds validity and can be utilized to verify our methodology.In the Angrist dataset, a child must be six years old within the current year to enroll in school under the U.S. Compulsory Education Law. In the U.S., the school year typically starts in August, meaning a child turning six in December can still commence their education in the same year. Consequently, a child born in the fourth quarter, such as December, can start school before reaching six. Conversely, a child born in the first quarter, like January, must wait until the autumn term after their sixth birthday to begin school. U.S. law mandates students must be at least 16 years old to legally drop out of school. Therefore, students dropping out at 16 may have varying years of education based on their birth month. For instance, those born between 1920-1929 have average educational years of 11.39, 11.44, 11.55, and 11 .57 for each quarter, respectively. Parents, when deciding to have children, seldom consider such subtle differences in birth months. Thus, the month of a child's birth, independent of other factors affecting educational levels like intelligence, family background, and environment, can be seen as a random assignment. This inadvertently creates variations in education duration based on birth month -akin to a randomized controlled trial where children born in the fourth quarter represent the "experimental group" with longer education, while those in the first quarter are the "control group" with shorter education. Hence, the birth quarter serves as an instrumental variable in this context.In this experiment, we continued to utilize synthetic datasets a and b for a classification study. This time, the labels were processed for binary classification. Specifically, we computed the probabilities for data points being classified into category 1 by applying a sigmoid function to the y values in the datasets; otherwise, the labels were assigned to category 0. Due to the inherent randomness in Question: Does the paper discuss the limitations of the work performed by the authors? Answer: [Yes] Justification: Please see Section 6. Guidelines:• The answer NA means that the paper has no limitation while the answer No means that the paper has limitations, but those are not discussed in the paper. • The authors are encouraged to create a separate "Limitations" section in their paper.• The paper should point out any strong assumptions and how robust the results are to violations of these assumptions (e.g., independence assumptions, noiseless settings, model well-specification, asymptotic approximations only holding locally). The authors should reflect on how these assumptions might be violated in practice and what the implications would be. • The authors should reflect on the scope of the claims made, e.g., if the approach was only tested on a few datasets or with a few runs. In general, empirical results often depend on implicit assumptions, which should be articulated. • The authors should reflect on the factors that influence the performance of the approach.For example, a facial recognition algorithm may perform poorly when image resolution is low or images are taken in low lighting. Or a speech-to-text system might not be used reliably to provide closed captions for online lectures because it fails to handle technical jargon. • The authors should discuss the computational efficiency of the proposed algorithms and how they scale with dataset size. • If applicable, the authors should discuss possible limitations of their approach to address problems of privacy and fairness. • While the authors might fear that complete honesty about limitations might be used by reviewers as grounds for rejection, a worse outcome might be that reviewers discover limitations that aren't acknowledged in the paper. The authors should use their best judgment and recognize that individual actions in favor of transparency play an important role in developing norms that preserve the integrity of the community. Reviewers will be specifically instructed to not penalize honesty concerning limitations. 


Decision: The paper makes an important contribution in understanding the causal attribution that drives an outcome variable. Rather than standard attribution scores that explain a predictive model, the authors propose an IV-based method that explains the true outcome variable. The method is motivated through a theoretical analysis of the effect of unobserved confounder on attribution scores. Reviewers agreed that the technical contribution and insight is note-worthy and therefore I'm recommending Accept. However, they also have two suggestions for improving the clarity of the paper:
1) It will be helpful to make it clear in the abstract and introduction that the goal is to explain the true outcome variable (not a given ML model)
2) The limitations in Appendix A and the limitations discussed in the rebuttal period should be included in the main text.

Reviews: review: paper_summary: paper_summary: This paper addresses the problem of estimating causal effects with feature attribution, by applying SHAP and Integrated Gradients (IG) to two-stage models with instrumental variables. On synthetic datasets, the proposed methods IV-SHAP and IV-IG can better recover the ground-truth causal effects, compared to SHAP and IG applied to predictive models. On real-world datasets, IV-SHAP and IV-IG are more aligned with prior knowledge than SHAP and IG.
strength_weakness: strength_weakness: - It is unclear how the expectations with respect to $\hat{M}_{\phi}(\tilde{x} | \bar{x}_t, \psi_t)$ in Equations (4) and (6) are empirically computed in the experiments. The (approximate) computation can potentially be expensive.
- What is the impact on IV-SHAP and IV-IG when the three assumptions of instrumental variables are violated? The paper lacks a discussion on this.
- What is the impact on IV-SHAP and IV-IG when the input features are correlated? The paper lacks a discussion on this.
questions: questions: - Besides SHAP and IG, is it appropriate to apply other feature attribution methods to a two-stage model in order to estimate causal effects?
- For the real-world datasets, what are the baseline values for IV-IG, and how are features removed in IV-SHAP?

Minor comments:
- There are typos in lines 172 and 178: $\frac{1}{\mathcal{N}}$ should be replaced with $\frac{1}{|\mathcal{N}|}$.
- I recommend omitting the $D$ in $\mathbb{E}_D$ of Proposition 2 for consistency.
- The use of Hoeffding's inequality in the proof for Lemma 7 require some boundedness condition on $h(u_i)$, which is reasonable to assume but should be mentioned nonetheless.
limitations: limitations: The limitations of the proposed approach itself seem adequately addressed.
paper_summary: paper_summary: The paper addresses the challenge of unobservable confounders in feature attribution methods, which can lead to misinterpretations. The authors propose a novel approach called "data-faithful feature attribution," which trains models free of confounders using instrumental variables (IVs) to ensure that feature attributions are faithful to the data generation process.
strength_weakness: strength_weakness: 1. The two-stage training process and the use of advanced techniques like IVs may be complex to implement and require significant computational resources, potentially limiting the method's practical usability.
2. The theoretical derivations are based on the assumption that the influence of unobserved confounders on the target features is linear. This assumption may not hold in all real-world scenarios, limiting the method's applicability.
questions: questions: 1. How can one reliably identify suitable instrumental variables in practice? What criteria should be used to ensure their effectiveness?
2. How does the method perform with non-linear effects of unobserved confounders? Have you conducted any experiments or simulations to demonstrate its robustness in such scenarios?
3. Have you conducted sensitivity analyses to evaluate the robustness of the method under different conditions, such as varying degrees of distributional shifts or the quality of IVs?
4. How does the proposed method scale with larger datasets and higher-dimensional data? Are there any benchmarks or performance metrics available to demonstrate its scalability?
limitations: limitations: Please refer to the above questions.
paper_summary: paper_summary: The paper introduces the problem that unobservable confounders in data can mislead and negatively affect the quality of feature attribution explanations. Feature attribution methods do not consider unobservable confounders which poses the risk of misinterpreting the attribution scores found, which could have consequences, especially in high-risk application domains where users rely on feature attribution to calibrate their trust in the model. To address the issue, the authors propose to decouple input features affected by unobservable confounding factors and those not. The latter are then used to train a confounder-free model, for which they show that feature attribution methods SHAP and Integrated Gradients are faithful to the true data generation process and robust.
strength_weakness: strength_weakness: - The paper is quite inline-math-text heavy, which limits its readability. I believe the readability and therefore clarity can be improved by restructuring math-heavy paragraphs (room for this could be made by e.g. shortening the introduction).
- Some things were unclear, please see the questions.
- Minor comments:
     - Lines 146-154: The long equations are slightly difficult to follow when written inline, and will be clearer if using \equation or \align.
     - Lines 156-159: Thought experiment results could be shown in a table, for more clarity.
     - Line 191: it’s -> it is
     - Please increase the font size in Figures 2 and 3 for better readability.
questions: questions: - I did not fully understand how you identified the unobservable confounders in the real dataset use case. How could practitioners identify unobservable confounders in their models? 
- Feature attribution methods explain a model output $f(x)$ by assigning attribution scores in the input instance $x$. I am unsure if I may have misunderstood your approach, but by fitting a model with a different data (one time with all variables, and once without variables affected by unobserved confounders), are you not training and explaining a different model? If this is the case, I see your contribution rather providing a training procedure on only instrumental variables and without confounding variables, leading to data-faithful attribution. Currently, the paper's contribution often reads as if you are proposing a new attribution method (i.e. referring to your approaches as IV-SHAP and IV-IG). Could you please clarify? 
- By excluding variables affected by unobserved confounders in training, does the model's predictive performance (e.g. measured in accuracy) suffer or does the task become easier for the model? Is there a trade-off between data-faithful feature attribution and model performance?
- How similar are the resulting attribution scores of IV-SHAP and IV-IG? It has previously been found that local explanation methods like feature attribution often disagree because they only locally approximate the model (Krishna et al., 2022: The Disagreement Problem in Explainable Machine Learning). Does removing the confounding variables lead to less disagreeing feature attribution? 
- In feature attribution, we are often interested in the features with the largest positive or negative attribution scores as we want to understand which features a model relied on most for its output. In this case, the faithful ranking of the scores matters more than the score itself. I'm curious, in your experiments, did the ranking of scores change with the instrumental variable approach? How much does being faithful to the data matter in practice if not? 
- Out of curiosity: What do you think, how would this approach / does this approach transfer to high-dimensional problems like image or text classification?

I will happily adjust the score when I understand the paper better.
limitations: limitations: Limitations and broader impact are adequately discussed in the appendix.

Comments: title: Response to the official comment of Reviewer 6T6P
comment: Thank you very much for your valuable comments. We appreciate your emphasis on the assumptions of instrumental variables and correlations of input features, which we agree are crucial factors for the effectiveness of IV-SHAP and IV-IG. We will make sure to highlight these points in the final version.
Thank you for clarifying my questions. I appreciate the authors for explaining that the performance of IV-SHAP and IV-IG can depend on (i) whether the three assumptions of instrumental variables are violated and (ii) the correlations between input features. I don't consider these as weaknesses specific to IV-SHAP and IV-IG, but rather general problems in causal inference. I encourage the authors to highlight the importance of (i) and (ii) so readers have such awareness when using IV-SHAP and IV-IG. For example, ablation studies with synthetic datasets can help readers assess the impact of (i) and (ii) on the effectiveness of IV-SHAP and IV-IG. Overall, the authors have addressed my questions, so I maintain my original score of 6.
Thank you very much for your comments. We appreciate your insights and will certainly follow your suggestions to make the scope and problem clearer in the final version.
Yes, thanks for the confirmation! I now understand the paper better and see the proposed method as a promising direction to recover causal effects in data. I'll raise my score from 4->7.

I would like to ask the authors to please make the scope and problem clearer in the final version of the paper, in order to avoid confusion with explanation methods that aim at explaining model outputs (like my case).
Thank you so much for your insightful follow-up questions and comments.

* Q1: Yes, your understanding is correct.

* Q2: Yes, your understanding is correct. Our proposed approach can be used as a tool for analyzing phenomena, making it particularly suitable in medical and financial fields.

* Q6: Yes, data-faithful concept attribution is a suitable scenario for applying the approach of confounders, and we appreciate this idea. The example mentioned in the rebuttal can be effectively modeled by treating the concept of lighting conditions as an unobservable confounder.

    Besides, we also think that data-faithful feature attribution is meaningful in the presence of pixel-level confounders. In adversarial attacks, a well-known example is adding invisible pixel-level noise to a panda image, leading to the model misclassifying the image as a gibbon [1]. When we conduct feature attribution on the image, those noisy pixels appear important, though humans cannot see their significance. If we could train a model to eliminate the influence of such noise, the attribution would better match the judgment of a human on the panda image and the shape of the panda. We think that there might be similar noise confounders existing in image classification. If we remove these confounders, the feature attribution results of the model would be more data-faithful, which might imply better model robustness. Of course, the example of data-faithful attribution in the image we mentioned above has some differences from our paper's focus. Firstly, the method for removing pixel-level confounders in images might not use instrumental variables. Secondly, data-faithful feature attribution, as we mentioned here, aims to help humans judge the reliability of a model, as in image classification, where model performance is usually the main concern. The performance-oriented goal is different from the goal of our paper on understanding the relationship between input features and target feature values. We believe this is an interesting topic that deserves more detailed thought in terms of motivation and methodology. 

We hope our responses have clarified your questions and helped improve your opinion of our work.

Reference

[1] Ian J. Goodfellow, Jonathon Shlens, and Christian Szegedy. Explaining and Harnessing Adversarial Examples. ICLR 2015.
Thank you for the clarifications and additional analysis for Q5. There are only a few follow-up questions and comments:

Q1: Could you please confirm whether my understanding is correct? Unobservable confounders are not required to be found. By being able to observe their effects, it is possible to identify the instrumental variables and make use of them to train a model with data-faithful feature attributions. 

Q2: Could you please confirm whether my understanding is correct? In this work, the focus is not on providing a way to achieve data-faithful feature attribution explanations for the sake of explaining model outputs, but on *understanding* the relationship between data and target values. The model trained on instrumental variables explains this relationship as feature attribution applied to this model are data-faithful. Can it be seen as an analysis tool for phenomena, e.g. a disease as in the example you mentioned? 

Q6: Just a comment: Feature attribution in image classification usually treats each pixel as a feature, but I like your connection to lighting conditions. Perhaps data-faithful concept attribution fits the approach of confounders as well, what do you think?