Explaining the Uncertain: Stochastic Shapley Values for Gaussian Process Models
Authors: Siu Lun Chau
Krikamol Muandet
Dino Sejdinovic
Venue: NeurIPS 2023 spotlight (2023)
Keywords: Gaussian Processes
Shapley values
Uncertainty Modelling
TL;DR: A new SHAP algorithm design for GPs that take into account the analytical covariance and result in explanations as Gaussian processes. A Shapley kernel is also proposed for predicting Shapley values of new observations.

Abstract: We present a novel approach for explaining Gaussian processes (GPs) that can utilize the full analytical covariance structure present in GPs. Our method is based on the popular solution concept of Shapley values extended to stochastic cooperative games, resulting in explanations that are random variables. The GP explanations generated using our approach satisfy similar favorable axioms to standard Shapley values and possess a tractable covariance function across features and data observations. This covariance allows for quantifying explanation uncertainties and studying the statistical dependencies between explanations. We further extend our framework to the problem of predictive explanation, and propose a Shapley prior over the explanation function to predict Shapley values for new data based on previously computed ones. Our extensive illustrations demonstrate the effectiveness of the proposed approach.

Explaining the Uncertain: Stochastic Shapley Values for Gaussian Process Models
Siu Lun Chau CISPA Helmholtz Center for Information Security 66123 Saarbrücken, Germany siu-lun. chau@cispa.de
Krikamol Muandet" CISPA Helmholtz Center for Information Security 66123 23 Saarbrücken, Germany
muandet@cissa.dde
Dino Sejdinovic*
School of Computer and Mathematical Sciences University of Adelaide Adelaide, South Australia, 5005 Australia
dino . www.einoiiicddelllieeee.. edu, au
Abstract
We present novel approach for explaining Gaussian processes (GPs) that can uti- lize the full analytical covariance structure present GPs. Our method is based on the popular solution concept of Shapley values extended to stochastic cooperative games, resulting explanations that are random variables. The GP explanations generated using our approach satisfy similar favorable axioms to standard Shap- ley values and possess a tractable covariance function across features and data observations. This covariance allows for quantifying explanation uncertainties and studying the statistical dependencies between explanations. We further extend our framework to the problem of predictive explanation, and propose Shapley prior over the explanation function to predict Shapley values for new data based previously computed ones. Our illustrations demonstrate the effectiveness of the proposed approach.
Introduction
Shapley values [T], a solution concept derived from cooperative game theory, are an essential tool in explainable AI that helps understand and interpret the prediction of complex machine learning models. In particular, SHapley-Additive Additive-e (SHAP) [2 algorithms have been applied in various fields, including finance B. healthcare 4 and robotics [5). SHAP algorithms compute Shapley values by formulating cooperative game payoffs based on the function obtained from the learning algorithm, e.g., regressor. This function is often treated fixed, resulting deterministic game payoffs. In contrast, we argue that uncertainty in model predictions also plays an equally important role for trustworthy machine learning models as it enables users to make more informed decisions (6), assess the model's confidence 77. and identify areas where more data or improved features are needed [8]. Since uncertainty in model prediction can propagate to downstream expla- nation, explainability tools should also help users better understand the model's confidence or lack thereof in its predictions, allowing them calibrate their trust in the explanations. By providing transparency into both the model's 's predictions and uncertainty. Al systems can become more reliable, interpretable, and trustworthy, facilitating broader adoption across domains where trust is imperative.
"Equal contribution.
37th Conference on Neural Information Processing Systems (NcurIPS 2023).

---

Gaussian processes (GPs) 99 are natural model class for developing an explanation method that can account for predictive uncertainty. While there exist model-speeific explanation methods for popular models such as LinearSHAP [10] for linear models, TreeSHAP [1]] for trees, DeepSHAP [ for deep networks, and AKHS-SHAL [2] for kernel methods, there a lack of Gi-specific Shapley value based explanation methods despite the ubiquity of GPs in machine learning. While feature importance can be studied via automatic relevance determination (ARD) lengthscales [|] of the covariance kernel, they only reveal global notions of importance which are not associated to the specific predictions, and are limited to the specific classes of kernels. Yoshikawa and Iwata [33 have proposed using GP with local linear regression component for interpretability, but this approach results in a specific interpretable GP model rather than a general GP explanation algorithm. Recently, Marx et al. [14] proposed to sample multiple realisations of the GP function and to apply model-agnostic explanation algorithm to each sample, thus getting a distribution over explanations. While this in principle can give uncertainties in explanations, it is computationally expensive and not tailored for GPs. It does not take advantage of the various important properties that GPs enjoy, such as the fact that posterior means are reproducing kernel Hilbert space (RKHS) [15] functions and that we have fully tractable covariance structure. These properties, as we later demonstrate, can lead to more effective estimation of Shapley values and help to understand the statistical dependencies across these explanations. Additionally, this covariance structure allows us to formulate explanations themselves as GPs and thus are able to predict explanations of unseen data.
To address this gap, we present the GP-SHAP algorithm. In Section 2 we contextualize our work by introducing the concepts of stochastic cooperative games and stochastic Shapley values. Both are studied thoroughly by game theorists [66 17). but have not been utilized by the explanation community. Section 3 characterizes the stochastic cooperative game induced by GPs through the of the conditional mean process [18 8. a GP of conditional expectations. We use the weighted least squares formulation of Shapley values [19 2 express explanations as multivariate Gaussians with tractable covariance structure. We then introduce the GP-SHAP algorithm as a method for estimating these stochastic Shapley values. In addition, we extend the GP-SHAP algorithm by combining it with the Bayesian weighted least squares approach by Slack et al. 20 The extension we call BayesGP-SSh SHAP- integrates two sources of uncertainty: GP predictive posterior uncertainty as well as the Shapley value estimation uncertainty (in the cases where Shapley values are too expensive to compute exactly). Section 4 expands beyond GP explanations, demonstrating the applicability of our framework to other predictive models. We consider the setting of predictive explanations and propose a Shapley prior to represent explanations for more general predictive models as GPs. This allows us to predict explanations for new data without rel relying on the standard Shapley value procedure and provide predictive uncertainty to more general explanations beyond GP models, such as explanations from TreeSHAP, DeepSHAP, and KerneISHAP, Illustrations of the effectiveness of GP-SHAP, BayesGP-SHAP, and the Shapley prior for predictive explanations are provided in All proofs in this paper are given in the appendix.
Section 5 and conclude the paper in Section 6
2 Stochastic cooperative games and their Shapley values
In this section, we review the concepts of stochastic cooperative games, denoted as s-games, and their corresponding stochastic Shapley values, referred to SSVs. These concepts provide the necessary language to introduce GP-SHAP in the coming section. It is important to note that game theorists have studied s-games in various forms [16] 17 21, but they have not been adequately introduced in the explanation literature. While Covert and Lee (2 briefly discussed s-games, they focused on computing deterministic Shapley values (DSVs) of the mean of the s-game, overlooking the inherent uncertainty in s-games and their potential application to uncertainty-aware explanations.
Problem formulation. The results presented below are analogous to the original work of Shap- ley's 's [I. but adapted to games with random outcomes. Formally, let l denote the player set and let 2 L2 (R) be a s-game with random variable payoff. We restrict our attention to stochastic payoffs with finite second moments as we aim to characterise the variance of the corresponding stochastic Shapley values. However, this is not necessary condition to show existence and unique- ness of the SSVs. We further require v(0) = dep where ón, is the Dirac measure at some constant to R. We introduce the following concepts to formalise the axioms for SSVs:
Definition 1 (Carrier). A carrier of v is any N C Ω such that v(S) v(N S)) for all S S Ω.
2

---

Definition 2 (Permutation s-game). Denote [1) the set of permutations on Ω. For E [I),, the induced permutation s-game is 1/ (=S) := ( for all S S 9, where TS the image of S under n. Definition 3 (Stochastic value allocation). Denote G the space of s-games for Qt We say Φ g + L2 (R) [E2] a value allocation of v that assigns to each player i E Q stochastic value Φs (v)).
The natural extension of the original axioms by Shapley are extended to s-games as follows:
1. (s-symmetry) For any II(S2), Pni @.(v).
2. (s-efficiency) For each carrier N of v, EIEN 4(v)= = =v(N).
3. (s-linearity) For any v.w E 9, @lv + w) o(v) + o(w) where addition of s-games are defined as (v+ (S) v(S) (S) for all S _ Ω.
Note that all equality happens at the random variable level. Unsurprisingly, there is a unique stochastic value allocation that satisfies these three axioms:
Theorem 4 (Stochastic Shapley values). The only stochastic value allocation & of v satisfying s-symmetry, s-efficiency, and s-linearity takes the following form,
(1(1)=
Σ cisl (v(SUi) - v(SS)
SCNMO
where N is the smallest carrier set of 91s = N
[s]
and (v((v) is the 4th SSV of s-game v.
Note that Equation (T) is equivalent to the DSVs, except it is written in terms of random variables. We emphasise that this result has been proven in Ma et al. [16] using top-down approach, i.e. starting with Equation ( and verify it satisfies the stochastic axioms and uniqueness. In the appendix, we offer a contrasting perspective where we mirror Shapley's original bottom-up derivation, i.c. began with the stochastic axioms and subsequently determined the unique solution.
2.1 Variances of stochastic Shapley values are not Shapley values of variance games
As ① is now defined by summing over weighted differences between random variables. we can analyse the corresponding mean and variances across the SSVs. In the following. denote 6: 9 + R as the deterministic Shapley value allocation where 9 {v s -> R} is the space of deterministic cooperative games, referred to as d-games from now on.
Proposition 5. Given the player set St let v be stochastic game, & a stochastic Shapley value allocation, and 6 deterministic Shapley value allocation. Suppose that [E[v] and V[v) are the corresponding mean and variance d-games, respectively. Then, E[o(v)] = e(B[p]), but V[o(v)] + ((V(v)) In particular, the SSV variance is given by
V(o:(1)] Σ EE
SENMI 9s19st (C(ysui, us Clvs ]s - Clvs, usrui + Clus,
where us v(S) and C is the covariance function between the stochastic payoff's.
Proposition 5 highlights the difference between variances of stochastic Shapley values, which capture the propagated uncertainties through the s- s-game, and deterministic Shapley values of variance games 23 24 i.e. deterministic game D(S) V[)(S) for all S S for some stochastic game v. Therefore, in the context of model explanations, variance-based based allocation techniques such as the Shapley effects (25 and the work of Fryer et al. [26). where they computed Shapley values on games constructed using the variance and coefficient of determination R2 respectively. cannot be used to capture uncertainty of explanations, as those approaches are themselves explaining the model variance instead. Besides stochastic Shapley values, it is possible to leverage the framework of coalition interval games [27] and interval Shapley values (28) to provide confidence intervals for feature explanations, as Napolitano et al. [29] demonstrated. However, their approach is based on confidence interval of predictive models, thus this frequentist approach is not applicable to Bayesian models such as GPs and the obtained uncertainty around explanations have very different interpretations as well.
At first glance, computing the variance of each stochastic Shapley value seems cumbersome as requires summing over all possible coalitions twice. However, there exists a compact expression for the full covariance matrix of stochastic Shapley values across players for stochastic games constructed using a Gaussian process model, which we demonstrate in the following section.
3

---

3 Explaining GPs with GP-SHAP
This section is dedicated to introducing GP-SHAP, where the key idea about efficiently computing stochastic Shapley values for stochastic games that models conditional expectations of Gaussian processes. Specifically, in Section 3.1 we review the process of constructing Shapley values for model explanations using deterministic cooperative games. We then proceed to formalize the stochastic cooperative game induced by pushing the GP functions through the d d-game and demonstrate that the resulting stochastic Shapley values follow a multivariate Gaussian distribution. Section (3. is dedicated to the estimation procedure and the introduction of two algorithms for computing GP explanations: GP-SHAP and BayesGP-SHAAl GP-SHAP provides explanations that incorporate the GP predictive uncertainty, while Bayes -HSP extends this further by accounting for the additional estimation uncertainty arising in the weighted least squares procedure. This extension is similar to how BayesSHAP introduced by Slack et al. 20 extends SHAP.
3.1 Formulating stochastic Shapley values for GP models
Let Y be random variables taking values in the d-dimensional instance space X C R and label space Y (could be in E or discrete) respectively with the joint distribution p(X, We use [d] := {I,....d d} to denote feature index set and S C [d] as the feature index subset. For supervised learning with GPs, the usual goal is to model P(Y X x) (x())) [9] where f is the GP function of interest and e is a problem-specife transformation, e.g., 1 is the identity map for regression and sigmoid transformation for classification problem. We then posit prior f f GP(O, k) where & X x X' + R is a covariance kernel, and compute the posterior distribution (or its approximation) p(f D), which is typically also a GP with posterior mean mn and kernel k.
Explaining deterministic with d-game. To explain deterministic function f. we create d-game by treating each feature as a player and constructing the game with f. The resulting DSVs then provide feature attributions that satisfy various favourable properties, such the sum of the attribution equals to the prediction made at that specific point, identical features receive the same attribution, and null feature receives zero attribution [2]. The d-game for local attribution an observation x typically involves computing the expected value of f(X) using a reference distribution xs) for different feature subsets S C [d].
r(X X5
( S) == [x() X5 - xs).
(2)
This approach, also known as removal-based explanation [30] has been used in various explanation algorithms 22 [1 These d-games determine the "worth" of features in by looking at the expectation after "removing" the contribution of other features in [d]\S through integration. The choice of reference distribution leads to explanations with different properties, such as improved locality of estimation (32), promotion of fairness (33, or the incorporation of causal knowledge [34] 351. In this paper, we focus on the scenario where r corresponds to the data distribution, i.c., r(X Xs) p(X Xs), it is one of the most frequently used approaches in the literature [36 37 12 31)
Although different reference functions result in different estimation procedures, our formulation of s-games and SSVs for GPs still holds.
Induced s-game by stochastic f. Now suppose f is a random function with a posterior distribution p( D). The corresponding stochastic cooperative game with players [d] can be defined analogously to the deterministic case (p(/ID) X x 2(d] + (TR) such that, for a given observation E X and feature subset S C [d], the stochastic payoff is
"P(rID) (x, S) Ex (S(X Xs = xs]
For brevity, we write (P()DD) as Vyy, bearing in mind that V/y is a random function since f is a random function. In partieular, when f is a GP with mean function in and covariance function k the stochastic payoff y is also a GP indexed by x and S.
Proposition 6 (Stochastic game yj as induced GP). Lei ~~ GP(n,, k) with integrable sample paths, i.e. fx fddxx < oo almost surely. The stochastic payoff function V induced by f is a Gaussian process with the following mean and covariance functions:
(x, 5) == (mx()) Xs ss= =xs].
hy ((x, (x', = Ex.xx [(X,X') Xg xs, X's =x's
(5)

---

In other words, the stochastic game V/ now assigns payoff distribution to each feature subset S, where the variability arise from taking conditional expectation on a random function f. We note that this GP is also known as the conditional mean process, previously studied in Chau et al. [8] 18]
Stochastic Shapley values as multivariate Gaussians for GPs. Given the s-game Vf, we can characterise the corresponding SSVs using the weighted least squares formulation of DSVs 19 2 For each coalition S [d], let 2j E {0, 1 be the binary vector such that 2, = if S and Z E {0, 10 xd the concatenation of all 24 vectors. When the context is clear, we use Z and [d] interchangeably. Let W R=xx2
x2" be the diagonal matrix with entries Wi = w((,)
d-1
for all subsets with size 0 < |Sy] d. When Sj has elements, the weights
(15y ]185 (188((a-153)) 18;1)
are set to to enforce the efficiency axiom.
Theorem 7 (Stochastic Shapley values of vz). Let V, be an induced stochastic game from the GP f (m, k) and denote vx (x, Vy (x, S2a)] the vector of stochastic payoff's across all coalitions, then the corresponding stochastic Shapley values Φ(vg(x, )) follows a d-dimensional multivariate Gaussian distribution,
o(vg(xx-)) (183))~~((AAN (N (.) with A Z WZ)-'Z'V W.
where E[vy] R²2 and v[v=] x2* of the payoffs.
are the corresponding mean vector and covariance matrix
The derivation is straightforward using the fact that the stochastic payoff's Vx are multivariate Gaussian random variables and the matrix A is obtained from the weighted regression formulation of Shapley
values. The resulting variance of ( (x) can then be interpreted as the uncertainties around the explanations, propagated from the posterior variance of f.
3.2 Estimation algorithms: the GP-SHAP and the BayesGP-SHAP
In this section, we introduce our main algorithm GP-SHAP and its variant BayesGP-SHAP to estimate the stochastic Shapley values for GP posterior distributions.
Formulating GP-SHAP. Given a set of observations D y: yX with GP prior T~GP(O, k)), we can obtain the posterior GP f D GP k) by using the Gaussian conditioning rule for regression or e.g. the variational Laplace approximation for classification. In fact, the conditional expectations of mean and covariance function from posterior GP can be estimated using conditional mean embeddings (38) without the need of explicit density estimations, following derivations of Chau et al. [8 18]. We provide technical discussion on conditional mean embeddings in the appendix.
Proposition 8 (Estimating vz). Given D = (X, and the posterior GP f D g k), the mean and covariance function of the stochastic cooperative game 1/y can be estimated as,
(x S) = b(x, mX)),
((x, (x', 8²)) = b(xx Kxxb(x',S"
(7)
where b(x, S) = (Kx.X + - (X Xs). (X) [m(x]),... "m(xx))) and ks X's x Xss + R is the kernel defined on the sub-feature space of X and we write (ss, Xs)= == (xs( X15)...,kssl respectively. The parameter X > 0 is a fixed hyperparameter to stabilise the inversion.
Xns)) and K Xxx and Kxx as the gram matrix of X using kernel k and k
It is worth noting that the estimation my coincides with the one deployed in RKHS-SHAY 12 since we are also utilising conditional mean embeddings for the estimation of conditional expectation of RKHS functions. addition, we obtain the analytical covariance function that used to estimate the covariance of the stochastic Shapley values in the following proposition:
Proposition 9 (GP-SHAP). Let the matrix A be defined as in Theorem 7 for the multivariate stochastic Shapley values can be estimated as,
The mean and covariance
$ (x, :) N (AB(x,
[d] AB(xx(]( Kxx (d)AT)
(8)
where B(x, [d]) = [b(x, (d]1).... b(x, [d]
5

---

The complete algorithm, along with a discussion of computational techniques deployed to reduce computation time, such as vectorisation across using tensor operations and incorporation of the sparse GP formulation to speed up computations of quantities from Propositions [|and g is provided in the appendix. is worth noting that subsampling coalitions to reduce computational cost while estimating Shapley values is also a standard approach in SHAP algorithm implementations (2]. Empirically, this procedure has been shown to give an unbiased estimate of the true DSVs, with the variance decreasing at rate of O(t), where É is the number of coalition samples [22]. However, this approach incurs an additional source of uncertainty, due to estimation, and we propose to incorporate that well into our framework by utilising the BayesSHAP approach proposed by Slack et al. [20].
Formulating BayesGP-SHAL To capture the estimation uncertainty, Slack et [20] proposed BayesSHAP and reformulated the WLS procedure as hierarchical Bayesian weighted least square and studied the corresponding posterior distribution over the deterministic Shapley values. We provide the hierarchical data generation process below, with an abuse of notations, for a deterministic f, we have
1 z, S, Φ 2 01 ~ ~(0,001)
e el >~ ~\ wz)(z)) o² ~~ ~
(9)
o
(10)
where w and are the weight function and binary vector respectively introduced in Theorem 7 and lo. o are two hyperparameters, typically set to small values to keep priors uninformative. Slack al. [20] showed that the posterior distribution on σ² and e follows scaled Inv- x² and normal respectively, due to their corresponding conjugacies with the likelihood.
Proposition (BayesSHAP 20
Given the data generation above, the posterior distribution on 0
and o² follows:
0² Ze, J. x, N(ANx (Z) (ZZ))))
(11)
0² Ze.fxx ~ Scaled-Inv- x² (a (6o+e, lood + (Vx)
(12)
lo +
where E is the number of coalitions S = - we sample uniformly from 2|4], Ze is the binary matrix representing S, and We is the corresponding weight matrix, and A = (Z W ZeZe)-LZ, We is the WLS matrix, Vx [D S1), ..., (x Se)] is the vector of deterministic payoffs, and
? [(Vx - ZeAov)) T - ZeA (A.Vx) (AAv))
(13)
measures the average weighted error in the regression and the norm of the mean explanations.
is important to highlight that while both GP-SHAP and BayesSHAP have notion of the "posterior of Shapley values", the two sources of uncertainty which these approaches capture are very different: GP-SHAP corresponds to the predictive uncertainty induced by the GP posterior p(ʃ D), whereas BayesSHAP is the uncertainty due to having to estimate Shapley values when their exact computation is infeasible due to the exponential number coalitions - and is in Slack et al. (20) proposed for a deterministic f. Nonetheless, by integrating the BayesSHAP posterior of Shapley values through the posterior GP p(f D), we can in fact incorporate both sources of uncertainty, leading to our second algorithm, BayesGP-SHAP. SHAP.
Proposition 11 (BayesGP-SHAP, Continuing from Propositions D aan 10 the posterior distribution of the stochastic Shapley values can be estimated using the Bayesian WLS approach as,
Φ 02 Z,, (A, (X),, A,X),,,B((x)) S) (x + We,)
where o² is sampled o2 Ze ~ Scaled-Inv- x² (lo + (n0 bott
We note that in the above proposition, instead of integrating p(o² Ze, f. x, D) with respect to the posterior GP, which leads to a complex scaled mixture of Gaussians, simplify the expression and construct a scaled inverse chi-square distribution with s² (E[v,,)) instead, which represents the error of the weighted regression with respect to the mean payoff's E(vx...
Conditionally on a², the posterior variance in BayesGP-SHAP is therefore the sum of the variance from GP-SHAP and BayesSHAP due to Gaussian conjugacies.

---

4 Predictive explanations using the Shapley prior
In this section, we move beyond standard GP explanations by formulating explanation functions of a broader class of models as Gaussian processes. We introduce Shapley prior over the space of vector-valued explanation functions : X Fd, which correspond Shapley values for arbitrary functions f. By treating previously obtained Shapley values as regression labels, we can predict explanations for new data without relying on the standard procedure to compute Shapley values. Our framework is not limited to explanations generated from SP-SH but can also be applied to other explanation methods such as TreeSHAP, DeePSHAP, or model-agnostic KernelSHAP. The proposed approach learns the direct mapping from X' to the space of Shapley values, without the need to access the underlying model f during training, which different from previous predictive approaches such as FastSHAP [39]. In contrast to prior work, Hill et al. [40] also considered fitting a GP regression model directly to explanation functions, but their focus is on incorporating the corresponding predictive uncertainty to explanations from black-box classifiers rather than predicting Shapley values for unseen data.
To achieve this, we leverage the key ideas behind the GP-SHAP algorithm, which uses the facts that GPs remain tractable under conditional expectations (to build the g-game) and linear combinations (to GP prior GP(O, k) instead
compute the stochastic Shapley values). By applying these properties of a posterior, we have an induced prior over the space of Shapley functions Φ.
Proposition 12 (The Shapley prior over Φ). The prior GP(O, k) and the corresponding stochastic game v (x. S) = [Xf(X) Xg = xs] induce a vector-valued GP prior over the explanation functions & ~ GP x.) where K : X x X + pdxd is the matrix-valued covariance kernel
n(x, x') = A(x) A(x'),
A(x) = (
(14)
where (x) [B[ X) Xs, x..). .......... B[M(, XSya = xs and the T sign refers to taking inner products in the RKHS of k.
By utilizing the Shapley function prior, we can apply our approach to predict explanations for a wide range of models, including trees, deep neural networks, or RKHS functions, by treating their explanations as noisy samples from this prior and using them as regression labels. Our approach is based on the perspective that there exists true data generating mechanism f. X + y and that any model f we use is an approximation of fr. Therefore, any explanations derived from f can be seen as an attempt to reveal the true explanations under f.. This perspective has also been adopted in the work of Chen et al. (36] and Marx et al. [14.. We present the predictive posterior below.
Proposition 13 (Predictive explanations as multi-outpu GPs). Given D' = {(x,; o.) (X, 2x) where Φi Id are the Shapley values computed under predictive model f and 4x = (oi..... on)" the predictive explanations for new data x' distributed as,
6(x') D& ~ N (ma(x'), x(x,xx x(x', X (x,x
(15)
where mx(x)= bx X)T veo(Oxx) bx XX = (Xxx + o3 - xx,, XxX), x'), Xxx is the gram matrix for kernel K of size nd x nd, x(x', X) [<(x', xI). xx(x' xn)] is of size d x nd and oz is the noise parameter for regression.
In fact, we see that the posterior mean (x') from 15) are Shapley values of the following payoff vector Vx,, computed based on the observed explanations Tx.
Proposition 14 (Posterior mean as Shapley values for payoff vector Vx.. The posterior mean
m (x)) corresponds to Shapley values for the payoff vector Vx, i.e., mx(x') = AAx', where
Vx D" V (xx) T (x.)) and a; E Id is the li,..., + (d - 1)] subvector of (Zxx
a () (x)(
5 Illustraions
We demonstrate the proposed approaches through three sets of illustrations. We first conduct an ablation study to examine how predictive and estimation uncertainty captured by our explanation algorithms vary under different configurations. We then move on to discuss various exploratory tools that can assist practitioners in utilizing the stochastic explanations. Finally, demonstrate the
7

---

bayesSHAP
GPSHAP
BayesGPSHAP
25% data, 50% coalition used
25% data, 100% coalition used
SSWS
100% data, 50% coaltion used
100% data, 100% coalition used
SSVs
Madinc
HouseAge
AveRaoms
AveBldrms
Population
Aveoccup
Latitude
Longitude
Medinc
HouseAge AveRooms
AveBodrms
Population
Aveoccup
Lathude
Langitude
Figure Ablation study on different uncertainties captured by GP-SHAP, BayesSHAP, and BayesGP- SHAP when computing local explanations (SSVs) using the California housing dataset [41]. 95% credible intervals around explanations are shown.
effectiveness of the Shapley prior for the predictive explanations problem. Our code is included the supplementary material, and we provide implementation details in the appendix. For all experiments, we pick the radial basis function (RBF) as our covariance kernel l for the GPs.
5.1 Ablation study on different notions of uncertainties captured
We conducted a comparison between GP-SHAP, BayesSHAP, and BayesGP-SHAP to demonstrate the differences between model predictive uncertainty and estimation uncertainty captured in the stochastic Shapley values. BayesSHAP is applied to the predictive posterior mean of GP to obtain its explanation. For this purpose, we used the California housing dataset [41] from the StatLLib repository, which includes 20640 instances and 8 numerical features, with the goal of predicting the median house value for California districts, expressed in hundreds of thousands of dollars. We trained our GP model using 25% and 100% of the data and calculated local stochastic explanations from GP-SHAP, Haye and BayesGP-SHAI using 50% and 100% of coalitions. The results, shown in Figure I demonstrate that the magnitude of BayesSHAP uncertainties (green bars) are uniform across features as it is designed to capture the overall estimation performance and not feature-specific uncertainty. In contrast, GP-SHAP (red bars) exhibits varying uncertainties across features due to the propagation of variation from the posterior GP f to the attributions. This allows practitioners to make more granular statements about the uncertainty around specific feature explanations. We also observe that increasing the number of coalitions and training data generally leads to a decrease in both estimation and predictive uncertainties, but the estimation uncertainty drops more significantly. BayesGP-SHH PHAP (yellow bars) consistently provides more conservative uncertainty estimation by considering both predictive and estimation uncertainties.
5.2 Exploratory analysis of the stochastic explanations
We propose several exploratory analysis methods to aid practitioners in comprehending the stochastic explanations generated for their downstream tasks. To this end, we use BayesGP-SHi -HAP to explain Gaussian process model trained on the breast cancer dataset (42). The dataset consists of 569 patients, 30 numeric features, and the objective is to predict whether a tumor is malignant or benign.
Local explanation. To visualize local explanations, we can plot the mean and standard deviations of stochastic Shapley values in a bar plot. This approach allows us to not only understand the degree of contribution a feature has to the prediction but also the corresponding credible interval as measure of explanation uncertainty. perimeter" features have similar contribution in terms of their absolute mean stochastic Shapley values. However, the model is much more uncertain about the former, allowing the user to calibrate their trust in model explanations.
(alure2addsplayy that both "worst fractal dimension" and "worst

---

Bayasc PSHAP
Top 10 explanations BayeSGPSHAP for patient
wenst forri simess on
worst fraccal dimordlon
pendeer
warst tocue
wond concave powns
wnor
testure
dimens
wonst permeter
concauity em
compactness
mean contave points
pertieter
Sidt of sceskine sove
absokus 55%s
oslis cner
conconity
Alous sie moon b0hs
0.6
0.8
stochastic
acpiamation
(a) Local explanation
Global explanation
Graphical model for top 15 explanations for patient
Corrlation between top 10 explanations for patient
Fyoncatty smo.
0.75
Minit mrngonae
prempcchess praa
pirucaa dineroonn arod
MACre Toygal dinwneer)
Pouuuue aryy
-132
Harrs symmenyy
-4.5
Sparsity
einels emmy
(c) Local explanation correlation
(d) Local explanation graphical model
Figure 2: Illustraiions of possible exploratory analyses utilising the BayesGP- SHAP covariance
structure, applied on the breast cancer dataset.
Heuristic global explanation as means of absolute stochastic Shapley values. It is common to compute the mean absolute DSVs as proxy global feature contributions 3 based computed local explanations. Although this heuristic does not result in any explanations that are themselves DSVs of any game at the global level, they provide quick summary to practitioners about the overall contributions. We can similarly compute the average of the absolute stochastic Shapley values, which are now distributed according to folded multivariate Gaussian distribution 44 with tractable covariance structure. However, it is important to highlight key distinction from the previous approach: when computing the mean of absolute SSVs, we consider the uncertainty in the local stochastic explanations, while directly calculating the absolute values of mean SSVs (equivalent to computing absolute DSVs) disregards this uncertainty. We believe the former approach is more appropriate than the latter as practitioners utilising Gaussian processes would want their global explanations to take the predictive uncertainty into account. The blue and yellow bars in Figure 2b depict these two approaches, respectively. Notably, observe that the blue bars suggest that the worst fractal dimension" feature is considered more influential than "worst concave points" because it accounts for the higher variability in the former feature's 's explanation. In contrast, the yellow bar overlooks this variability, leading to a different conclusion.
Correlations and dependencies across explanations. As we have access to the explanation covariance, we can explore their correlation and visualise them as Figure [2 For instance, the stochastic Shapley values of "concavity error" are less correlated with other features. Moreover, explanations are multivariate Gaussians, we can build an undirected Gaussian graphical model [45) using the precision matrix (inverse of the covariance matrix) to study the independence across the local stochastic explanations. We visualise the corresponding graphical model for explanations from patient 1 in Figure 2d by setting sparsity threshold of 90%, following the approach in [46).
5.3 Predictive explanations
Finally, we showcase the effectiveness of our Shapley prior by comparing the regression per- formances of multi-output GP using the Shapley prior with multi-output random forest and neural networks when predicting withheld explanations generated from GP-SHAP, Tree- SHAP, and DeepSHAP, respectively. We use the diabetes dataset [47] from the UCI repos- itory, which contains 442 patients with 10 numerical features and the goal is to predict dis- case progression. We first train GP. random forest, and a neural network, to obtain the

---

subsequent explanations from GP-SHAP, TreeSHAP, and DeepSHAP. Next, we feed 70% of the explanations to our predictive model as training data and the remaining 30% as test data.
We repeat this process over 10 seeds and compute the root mean squared error between the predicted explanation and the exact explanations for respective models, as shown in Figure 3 We observe that our GP-model with the Shapley prior consistently outperforms random forest and neural network model for pre- dicting all three types of explanations. This demonstrates that the inductive bias that comes from our choice of covariance structure и allows to build more accurate predictive explanation models. We also observe that the overall average prediction error for explanations generated from GP-SHAP is lower than that of TreeSHAP [1T] and DeepSHAP 2 suggesting that GP-SHAP produce explanations that are easier to learn. This follows our intuition GPs are typically smoother functions than trees and neural networks.
Secanot on Rend gortere
17.--
Shapley
Random forest
Neural Network
7.5
0.0
GP-SHAP
TreeSHAP DeepSHAP
Sources explanations
Figure 3: Predictive performance of using Shapley prior predict ex- planations generated from different explanation algorithms the dia- betes dataset.
6 Discussion
In this work, we presented a novel and principled approach for explaining Gaussian process models using stochastic Shapley values. The proposed algorithm GP-SHAP, and its variant BayesGP-SHAA allow practitioners to consider both predictive and estimation uncertainties when reasoning about the explanations and calibrate their trust in the model accordingly. Furthermore, we consider the setting of predictive explanations, where we introduced a Shapley prior over explanation functions, enabling us model and predict Shapley values for wider range of predictive models. We demonstrated the effectiveness of our methods through a number illustrations and discussed various exploratory tools for practitioners to analyse the stochastic explanations.
Limitations and future outlook. While the general framework of using stochastic Shapley values for stochastic explanations can be applied beyond Gaussian process models, the specific estimation algorithms presented in this paper are tailored to GPs and may not be directly applicable other probabilistic models like Bayesian neural networks. These alternative models would require different estimation procedures as the resulting explanations would not GPs anymore. Another promising avenue for future research involves leveraging the uncertainty obtained from predicting explanations and incorporating it into Bayesian optimization
guiding experimental design. For example, this could allow practitioners explore data regions that produce significant Shapley values for specific features. At last, we would like to highlight that although SHAP is popular, as an explanation tool contains several limitations, specifically Kumar et al. 48 have shown that SHAP violates several expected properties of explanations from human-centrie perspective.
Acknowledgement
The authors would like to thank Jean Francois Ton, Shahine Bouabid, Jake Fawkes, and Simon Föll for insightful discussions. The authors are also indebted to anonymous reviewers whose feedback has significantly improved the manuscript.
10

---

References
[1] Lloyd S Shapley. A value for n-person games. Contributions to the Theory of Games, 2(28): 307-317, 1953.
[2] Scott M Lundberg and Su-In Lec. A unified approach to interpreting model predictions. In Advances in neural information processing systems, pages
2017.
[3] Philippe Bracke, Anupam Datta, Carsten Jung, and Shayak Sen. Machine learning explainability in finance: an application to default risk analysis. Bank of England Working Paper, 2019.
[4] Konstantin D Pandl, Fabian Feiland, Scott Thiebes, and Ali Sunyaev. Trustworthy machine learning for health care: scalable data valuation with the shapley value. In Proceedings of the Conference Health, Inference, and Learning, pages 2021.
[5] Sindre Benjamin Remman, Inga Strümke, and Anastasios M Lekkas. Causal versus marginal shapley values for robotic lever manipulation controlled using deep reinforcement learning. In 2022 American Control Conference (ACC), pages 3683-2690 IEEE, 2022.
[6] Umang Bhatt, Javier Antorán, Yunfeng Zhang, Q Vera Liao, prasanna Sattigeri, Riccardo Fogliato, Gabrielle Melançon, Ranganath Krishnan, Jason Stanley, Omesh Tickoo, et al. Un- certainty as form of transparency: Measuring, communicating, and using uncertainty. In Proceedings of the 2021 AAAT/ACM Conference on AI, Ethics, and Society, pages 401-413, 2021.
[7] Eyke Hüllermeier and Willem Waegeman. Aleatoric and Epistemic Uncertainty in Machine Learning: An Introduction to Concepts and Methods. Machine Learning, 110(3):457-55 March 2021. ISSN 0885-6125, 1573-0565. doi: 10. 1e.11000ss1990 00077//10994-02--554446 URL h //arxiv. org/ abs/ 1910. 09457 arXiv:1910.0945 [cs, stat].
[8] Siu Lun Chau, Dean-Francois Ton, Javier González, Yee Teh, and Dino Sejdinovic. BaY Uncertainty Quantification for Causal Data Fusion. In Advances in Neu-
Curran Associates, ral Information Processing Systems, volume 34, pages pttps///ppppee/////////aaaaa pptp:///p////aslll Inc., 2021. URL http:/ttsss/prrceeeiinggi// mew.:/eeeedddinnennneeeee..eeeee cttpppp..ccpppppe
wwwpppppppppppeeeeetee files/ (ca2775aaS01aa1115aaalllaaaaaaa Absssaddddc----raacccccaaaaaaaa- html
[9] Carl Edward Rasmussen, Christopher KI William,, et al. Gaussian processes for machine learning, volume Springer, 2006.
[10] Erik Štrumbelj and Igor Kononenko. with feature contributions. Knowledge and information systems, 41(3):647-665, 2014.
Explaining prediction models and individual predictions
[11] Scott M Lundberg, Gabriel G Erion, and Su-In Lee. Consistent individualized feature attribution for tree ensembles. arXiv preprint arXiv:1800 1802.03888 2018.
[12] Siu Lun Chau, Robert Hu, Javier Gonzalez, and Dino Sejdinovic. RKHS- -SAP: Shapley values for kernel methods. Advances in Neural Information Processing Systems. 35:13050-13063, 2022.
[13] Yuya Yoshikawa and Tomoharu Iwata. Gaussian Process Regression with Local Explanation. IEEE Transactions on Neural Networks and Learning Systems, pages 1-15, 2021. ISSN 2162- 237X, 2162-2388. doi: 10. 109/TNNLS.2000 2.21...222.31332544 URL http://arxiv/ org /bs/
2007 .06699 arXiv:2007.011 11000001666 669 [cs, stat].
[14] Charles Marx, Youngsuk Park, Hilaf Hasson, Yuyang Wang. Stefano Ermon, and Luke Huan. But are you sure? uncertainty aware perspective explainable ai. In Francisco Ruiz, Jennifer Dy, and Jan-Willem van de Meent, editors, Proceedings of The 26th International Conference on Artificial Intelligence and Statistics, volume 206 of Proceedings of Machine Learning Research, pages 7375-7391. PMLR, 25-27 Apr 2023. WRl.hhttss//prrceeee mlr. wttp://vv0mmmaaaae html
[15] Vern Paulsen and Mrinal Raghupathi. An introduction the theory of reproducing kernel Hilbert spaces, volume 152. Cambridge university press, 2016.
11

---

[16] Ying Ma, Zuofeng Gao, Wei Li, Ning Jiang, Lei Guo, et al. The shapley value for stochastic cooperative game. Modern Applied Science, 2(4):1-76, 2008.
[17] Jinwu Gao, Xiangfeng Yang, and Di Liu. Uncertain Shapley value of coalitional game with application to supply chain alliance. Applied Soft Computing. 56:551 July 2017. ISSN 15684946. doi: 10. wswo//.sss............oooooooo www........06.66 retrieve/pie pttppp////51994111500000
URL
bew.akkiinguu...seiieerr... elsevier com/
[18] Siu Lun Chau, Shahine Bouabid, and Dino Sejdinovic. Deconditional Downscaling with Gaussian Processes, In Advances in Neural Information Processing Systems, volume 34, pages 17813-17825. Curran Associates, Inc., 2021. URL httpp: http:/s.nppps.cccpppppppeppp htw.:/9adddd-aattttthhaheaee
[19] A Charnes, B Golany, M Keane, and J Rousseau. Extremal principle solutions of games in characteristic function form: core, chebychev and shapley value generalizations. Econometrics of planning and efficiency, pages 123-133, 1988.
[20] Dylan Slack, Sophie Hilgard, Sameer Singh, and Himabindu Lakkaraju. Reliable Post hoc Explanations: Modeling Uncertainty in Explainability, November 2021. URL http://axxiv
org/abs/2000 05030
arXiv:2008.05030 [cs, stat].
[21] Boyang Dai, Xiangfeng Yang, and Xiaoyue Liu. Shapley Value of Uncertain Coali- tional Game based Hurwicz Criterion with Application to Water Resource Allocation. Group Decision and Negotiation, 31(1):241-263 February 2022. ISSN 0926-2644, 1572- 9907. doi: 10. 10111000//10722--2---0776--- $10726-001-09767--6
htRpsttts:://ne springer com 1007/
[22] Ian Covert and Su-In Lee. Improving KerneISHAP: Practical Shapley Value Estimation
Using Linear Regression. In Proceedings of The 24th International Conference Ar- rificial Intelligence and Statistics, pages 3457-3465. PMLR, March 2021. URL https ISSN: 2640-3498.
/proceedings.mmr. wtt.p//...peee//11/////eeeerrr/e html
[23] Riccardo Colini-Baldeschi, Marco Scarsini, and Stefano Vaccari. Variance Allocation and Shapley Value. Methodology and Computing in Applied Probability, 20(3)9919-933, September
2018. ISSN 1387-5841, 1573-7713, doi: 10.1007/s1 1009-016-994000 9540-5,. URL TRRLhttp://// http:://tatttt /arxiv org/abs/1606 09424
arXiv: 1606.0944444 [cs, math, stat].
[24] Marcello Galeotti and Giovanni Rabitti. On the comparison of Shapley values for variance and standard deviation games. Journal of Applied Probability, 58(3):609 620, September 2021.
ISSN 0021-9002, 1475-6072. doi: 10. 10.10/0//pp... 106. URL https: :/tp::///wwwww. cambridge
[25] Art B Owen. Sobol and shapley value. SIAM/ASA Journal on Uncertainty Quantification, 2(1):245-25) 2014.
[26] Daniel Fryer, Inga Strümke, and Hien Nguyen. Shapley value confidence intervals for attributing variance explained. Frontiers Applied Mathematics and Statistics, 6:587199, 2020.
[27] Rodica Branzei, Oana Branzei, S Zeynep Alparslan Gök, and Stef Tijs. Cooperative interval games: survey. Central European Journal of Operations Research, 18:397-411, 2010 2010.
[28] SZ Alparslan Gök, Rodica Branzei, and Stef Tijs. The interval shapley value: an axiomatization. Central European Journal of Operations Reseurch, 18:131-140, 2010.
[29] Davide Napolitano, Lorenzo Vaiani, and Luca Cagliero. Learning Confidence Intervals for Feature Importance: A Fast Shapley- based Approach. Workshop Proceedings of the EDBT/ICDT 2023 Joint Conference (March 28-March 31, 2023, loannina, Greece).
[30] Ian Covert, Scott Lundberg, and Su-In Lee. Explaining by removing: A unified framework for model explanation. Journal of Machine Learning Research, 22(209):1 2021.
[31] Robert Hu, Siu Lun Chau, Jaime Ferrando Huertas, and Dino Sejdinovic. Explaining preferences with shapley values. Advances in Neural Information Processing Systems, 35:27664-27677, 2022.
12

---

[32] Sahra Ghalebikesabi, Lucile Ter-Minassian, Karla DiazOrdaz, and Chris C Holmes. On locality of local explanation models. Advances in Neural Information Processing Systems, 34, 2021.
[33] Masayoshi Mase, Art B Owen, and Benjamin B Seiler. Cohort shapley value for algorithmic fairness. arXiv preprini 1rXiv:2105.07168, 2021.
[34] Christopher Frye, Ilya Feige, and Colin Rowat. Asymmetric shapley values: incorporating causal knowledge into model-agnostic explainability. arXiv preprint arXhv:1910.06338, 2019.
[35] Tom Heskes, Evi Sijben, Ioan Gabriel Bucur, and Tom Claassen. Causal shapley values:
Exploiting causal knowledge to explain individual predictions of complex models. Advances in neural information processing systems, 33:4778-4789, 2020.
[36] Hugh Chen, Joseph D Janizek, Scott Lundberg, and Su-In Lee. True to the model or true to the data? arXiv preprint arXiv:2006.1622: 2020.
[37] Christopher Frye, Damien de Mijolla, Laurence Cowton, Megan Stanley, and Ilya Feige. Shapley-basee explainability on the data manifold. arXiv preprint wrxiv:2006.01272. 2020.
[38] Krikamol Muandet, Kenji Fukumizu, Bharath Sriperumbudur, Bernhard Schölkopf, et al. Kernel mean embedding of distributions: A review and beyond. Foundations and Trendy@ in Machine Learning, 1011,)): 2017.
[39] Neil Jethani, Mukund Sudarshan, Ian Covert, Su-In Lee, and Rajesh Ranganath. FaPt Real-Time Shapley Value Estimation, March 2022. URL httt:://arxiv.org/bbell
07436
arXiv:2107.07338 [cs, stat].
[40] Davin Hill, Aria Masoomi, Sandesh Ghimire, Max Torop, and Jennifer Dy. Boundary-Aware Uncertainty for Feature Attribution Explainers, February 2023. URL http://arxiv..om org/ 0bs/2210.000 02419 arXiv:2210.02419 [cs].
[41] Kelley Pace and Ronald Barry. Sparse spatial autoregressions. Statistics Probability Letters, 33(3):291-297, 1997.
[42] W Nick Street, William Wolberg, and Olvi L. Mangasarian. Nuclear feature extraction for breast tumor diagnosis. In Biomedical image processing and biomedical visualization, volume 1905, pages 861-870 SPIE, 1993.
[43] Ian Covert, Scott Lundberg, and Su-In Lee. Understanding global feature contributions with additive importance measures, Advances in Neural Information Processing Systems, 33, 2020.
[44] Raymond Kan and Cesare Robott.. On moments of folded and truncated multivariate normal distributions. Journal of Computational and Graphical Statistics, 26(4):930-934, 2017.
[45] David Edwards. Introduction to graphical modelling. Springer Science & Business Media, 2012.
[46] Salar Fattahi and Somayeh Sojoudi. Graphical lasso and thresholding: Equivalence and closed-form solutions. Journal of Machine Learning Research, 20(10):1-44, 2019. URL
http://jmlr.. httg:papers/v2//111 1000001715511 501 html
[47] Bradley Efron, Trevor Hastie, Iain Johnstone, and Robert Tibshirani. Least An- gle Regression. The Annals of Statistics, 32(2), April 2004. ISSN 0090-5364. doi: 10. 00,12114/09053344000000000 arXiv.math/0400456.
URL. http://arre: org/ arggabs/mmaaaa math//406h/0666
[48] I Elizabeth Kumar, Suresh Venkatasubramaniin, Carlos Scheidegger, and Sorelle Friedler. Prob- lems with shapley-value value-based explanations feature importance measures. In International Conference on Machine Learning, pages 5491-5500. PMLR, 2020.
[49] Michalis Titsias and Neil D Lawrence. Bayesian gaussian process latent variable model. In Proceedings of the thirteenth international conference on artificial intelligence and statistics, pages 844-851, JMLR Workshop and Conference Proceedings, 2010.
13

---

[50] Le Song, Kenji Fukumizu, and Arthur Gretton. Kernel embeddings of conditional distributions: unified kemeel framework for nonparametric inference graphical models. IEEE Signal Processing Magazine, 30(4))98-111 111,,, 2013.
[51] Mauricio A. Alvarez, Lorenzo Rosasco, and Neil D. Lawrence, Kernels for Vector-Valued Func- tions: Review, April 2012. URL http://arxiv, hwtp:/orr/bbbbbbbabee/ee 6251 arXiv:1106.622. [es, math, stat].
[52] Bernhard Schölkopf, Ralf Herbrich, and Alex Smola. A generalized representer theorem. In Computational Learning Theory: 14th Annual Conference on Computational Learning Theory, COLT 2001 and 5th European Conference Compulational Learning Theory, EuroCOLT 2001 Amsterlam, The Netherlands, July 16-19, 200 Proceedings 14, pages 416-426. Springer, 2001.
14

---

SUPPLEMENTARY MATERIAL: Explaining the uncertain: Stochastic Shapley values for Gaussian process models
A The GP-SHAP algorithm and discussion on computation techniques
We present the complete algorithm for both GP-SHAP and BayesGP-SHAP in Algorithm
|
Algorithm 1 GP-SHAP/ BayesGP-SHAP
Input:
Posterior mean function ri, posterior covariance function k, inducing locations X, expla- nation instances X, number of coalition samples nz, hyperparameter X, Ro, o- base kernel k algorithm algo,
1: Compute nı number of inducing location, nuumerr of explanation instances, d = number of features.
2: Compute Cholesky decomposition on posterior covariance LL
K xx
3: Sample coalitions S {S1, .... Sn²] from [d]. build binary matrix Z {0, 1]nzxe from S, and compute weights W diag|w1,-.. ..... Waz] with w,
d-1
(1e 15.
4: Compute A Z WZ W
Shape: d x nz
5: Compute B X S) [K x R + (Ks Xs) for S in S]
Shape: az X 11
6: Compute Q where Qi,. 2, B(X, S)
Shape: f2 X n. X 11
7: Compute R where R,A,,. E, Ar, QiAi
> Shape: X T. X n
8: Compute V where V Vi,,,k.n Eje Ri,jA Rum.n
b Shape: x d
n
9: Compute E where Ei, E, B(X, () (),,,(())
Shape: 72
10: Compute 4 = AE
The mean stochastic Shapley values of shape d
11: if algo = GP-SHAP then
12: return mean explanations T and covariance V between d features and instances
13: else if algo = BayesGP-SHAA then
14:
Compute s² diag (EE Z W(E ZΦ) + diag( + TT)
Shape: X 1
Sample ²² from Scaled- Inv- x² no + nz. no to
nao +ngs
Shape: n.
16:
return mean explanations and covariance + (Z WZ) -12
17: end if
Computational considerations. In terms of computational complexity, one of the most demanding operations in the algorithm is the computation of conditional mean embeddings in step Instead of naively inverting an n X n matrix, which would have computational cost of O(n²), ), we employ the conjugate gradient method reduce the computation of the conditional mean embedding component to On²a), where a & n represents the number of conjugate gradient iterations. Additionally, to further reduce runtime,
utilize the variational sparse GP model 49 This model learns set of
inducing locations X with a size of n. & n, which can be reused for the estimation of conditional mean embeddings in the algorithm. Consequently, the computation of the conditional expectation is O(ffa). Another computational burden arises from the computation of the full covariance matrix across d features and n. instances, which requires storage of n²d matrix. However, since the full covariance matrix can be factorized into the component from step of the algorithm, we can store this low-rank component and compute covariances between specific instances when necessary. worth noting that this decomposition of the covariance matrix allows us avoid redundant computations when computing the covariance component, as we no longer need to iterate over all possible coalitions twice. Finally, we can further speed up our computational by parallelising computation across the sub-sampled coalitions in step 5.
reduced from (n)aa)
15

---

B Proofs and derivations
B.1 Section 2 proofs: Stochastic Shapley values
We include the full proof of the derivation stochastic Shapley values for completeness. The proof analogous the original work of Shapley's [ but extended to random variable payoffs. Ma et al. [16] has also proved the same theorem but used a different proving strategy. They started with the solution and showed satisfies the axioms and then prove uniqueness, whereas the following proof starts from the characterisation of s-games and derive the solution from a bottom-up fashion.
To facilitate the proof, first introduce the concept of stochastic symmetric game.
Proposition 15 (s-symmetric games). Let C be a real-valued random variables, then the symmetric game (5) == C1(R C S] gets a stochastic shapley value as,
(s(vc.))= C
(16)
where [R].
Proof. Take any i,j E R, pick a permutation æ E (TU so that R and wi = j. so the induced game TVC.R R. and therefore by the s-symmetry axiom,
dy(Vc,n) (VC.n))
(17)
Now by the s-efficiency axiom,
C=VC,R(R) E =v)( (VC() jen
(18)
for any i R.
Now we can characterise the form of any stochastic game as follows:
Proposition 16. All s-games with finite carrier can be written linear combination of s-symmetric games,
V
Σ
Ven(w),R
(19)
RCN,RØ0
where
Cnlv) = (-1)-=u((T)
(20)
TER
Proof.
We start by verifying
MSS) = Σ (eR(v),R((S)
(21)
RCN
holds for all S C U, and for any finite carrier N v. IfS N. then we can rewrite the expression
as,
v(S) Σ Σ (-1)"
(22)
RESTER
Σ Σ (-1)"-tv(I)
(23)
TESTERCS
Σ ((T) Σ (-1))
8 - r r-t)
(24)
TES
-(S)
(25)
where in the last equation we used the fact that (-)r--1)-((--1 is binomal expansion of (1 + (-1))=-1, therefore the only non-zero expression is when = 8.
16

---

We can now prove the uniqueness of stochastic Shapley values. Theorem 4 (Stochastic Shapley values). The only stochastic value allocation & of satisfying s-symmetry, s-efficiency, and s-linearity takes the following form,
01(1)
Σ
91s) (v(SUi) -
(I)
SCNVO
(1 and (v( is the tth SSV of s-game
where N is the smallest carrier of 91s)
[
[s]
Proof. First, let us denote
(s()) Σ (-1)-
RCN SU
Applying the s- s-Iinearity axiom on ø to the characterisation of y from the previous propositions leads us the following.
0i(v) = di
Σ
"Ca(v),
(26)
REN
=
E
(i(VCR(U))))
(27)
RC
Σ
en(y)=
(28)
RCN,ER
Σ
((11)-(((5)
(29)
REN
SER
Σ
Σ
(-1)-~u(S))
!
(30)
SCN RCN
SUGJER
Σ (((S)b(S)
(31)
SCN
Σ
((S))(S)+ + +7 - (i))a(S (()))
(32)
SCN
Σ
74(5) -
(33)
SCN
ES
(s 1)l(n - s)! (V(S) - (S (i)))
(34)
n!
SCN
Σ
Cis (U(SU))-((
(35)
SCNYO
where in 332 we used the following observation: given i $ 5' C N, and S = S" U (i}, then (s)) = -
satisfies uniqueness by construction.
Proposition 5. Given the player set St let v be stochastic game, a stochastic Shapley value allocation, and ó deterministic Shapley value allocation. Suppose that [E(v) and V[y) are the corresponding mean and variance d-games, respectively. Then, [[((v)] = ((E[p), but V[o(v)) + o(V(v)). In particular, the SSV variance is given by
V[o.(v)] = Σ Σ cssss (Clvsur, vsull Clysun us:
Clps, us + Clrs.
SENVE) S')
where us v(8) and C is the covariance function between the stochastic payoffs.
17

---

Proof. The equivalence between mean of stochastic Shapley values and deterministic Shapley values of mean game is trivial show leveraging the linearity of expectation. The variance of V ( can be shown by repeatedly applying the standard identity V[X + YJ V + V[Y] + 2C X YJ for random variables X, Y. Now consider the deterministic Shapley values of variance game V(u).
Σ SC
cis
(36)
Comparing to the expression of (z) from the lemma,
V Σ Σ
9s19s (Ciusue, usuul Ciusua, Clvs, usruil + Clss)
SCNMA S'CNM()
if assume mutual independence across all payoff random variables, leading to Civ(S U i), v(S) = 0 for all S, we still would not subtract but instead sum the variance ([v(S U i) and
[p(()). Therefore the variances of stochastic Shapley values is not the same the deterministic Shapley values of the variance game.
B.2 Section 3.1 proofs on the stochastic Shapley values for induced stochastic game from GP
Proposition (Stochastic game V as induced GP). Let f ~ GP le) with integrable sample paths, i.e. Ix fddy < oo almost surely. The stochastic payolf function Vj induced by is a Gaussian process with the the following mean and covariance functions:
(x = Ex (m) Xg = =xs],
(4)
hi (x, = Ex,xx
x(X,X') Xs Xs, X's xx's
(5)
Proof.
This is a direct application of Chau et al. 18 Proposition 3.2] to the distribution P(X
Xs
xs).
Theorem (Stochastic Shapley values of vy). Let Vj be an induced stochastic game from the GP I ~ Gm, k) and denote vx := (x, Si)... vy (x, S2d) T the vector of stochastic payoffs across all coalitions. then the corresponding stochastic Shapley values ((v,(x,, 1) follows d-dimensional multivariate Gaussian distribution,
Φlv,(x,)) ~ ((AE AVV.aA)) with A = Z WZ) T W.
where E[va] E R²2 and V[vz] E x22" are the corresponding mean vector and covariance matrix of the payoffs.
Proof. Recall from Lundberg and Lee 2 Theorem 2]. for deterministic Shapley values, given deterministic payoff vx for all 2d coalitions, the expression of Shapley values for each i E [d].
Σ Css
(37)
SCION()
can be written compactly as the following vector,
Hx AVx
(38)
We can therefore similarly write down the form of the stochastic Shapley values using this linear operator A, acting now on a vector of random variable output stochastic payoff vector Vx-
Avx
(39)
Nonetheless, as Proposition 8 implies that Vx is multivariate Gaussian, therefore ox is also multivariate Gaussian with mean and covariance the following,
Vx N AV
(40)
18

---

B.3 Section 3.2 proofs on estimation
To proceed, we first introduce the concepts of conditional mean embedding as a tool to estimate conditional expectation of functions living in their corresponding RKHSs,
Definition 17 (Conditional mean embedding [8)). Let X, be random variables and k X + X R kernel on X, then we define the following as the conditional mean embedding of p(X Y = y),
=x :- f (dPP))
(41)
Proposition 18 (Conditional Mean estimation). For random variable X.Y Y, and kernel k X + X - IR on X and kernel 1: y - 3 - R on V. Given observations D = {X,y), the empirical conditional mean embedding can be estimated as
PX|Y= = V( (Iyy + - +(X,),
(42)
where y) = v1).......... .,, yn)) and k(,, X) [k(,, xI)....... k(,, xx)] the parameter X > 0 is there to stablise the inversion. Now for Hi. the conditional expectation can then be estimated as,
B(S(X) = -
(43)
= -I(,,)((ssv + -f.
(44)
where f = (f(x)),....) (xx)]
Proof. This is standard result from literature, please read Song et al. 50 Muandet et al. 38 for more details.
Now we can apply these propositions to estimate the mean and covariance functions of the induced stochastic game from GP,
Proposition 8 (Estimating vy). Given D y) and the posterior GP f D GP k), the mean and covariance function of the stochastic cooperative game V/ can be estimated as,
(x,(x, S) S) mx(),
(xx, (x,, S')) b(x, S) K Xx'
(7)
where b(x, S) == (KxsXs + +
m(X) = [m(x)).... ,m,,))] and ks Ns Xss + is the kernel defined on the sub-feature space of X and we write (552x,, [xg(Xs: (is)..,.).....ee (xxs,,, xns)). and Kxx and Kxx as the gram matrix of X using kernel l and respectively. The purameter >> 0 is fixed hyperparameter to stabilise the inversion.
Proof. Without loss of generality, we will demonstrate this proposition with m. h obtained via standard GP regression, i.e.,
mxx) k(x, X) (Kxx 021)
(45)
M(x, x') = k(xx x') - X) + -
(46)
Starting with the mean function,
E(m(X) Xs xs] Ex( X) K xx + 021) ly Xg= x
(47)
= ())(,X)) + 0 1 y. xs=xs) H H2.
(48)
We can replace the population conditional mean embedding with the empirical version, and expand,
Elm(X) X xs] (A(,, X) + +0²1)= ;xixxxxxx No.
(49)
kXs s, X
K xx + AI) "xxx +011
(50)
=b(x, (X).
(51)
Analogously, the conditional expectation of the posterior covariance function, i.e., E[i(X,, Xg = Xs, X = xzy, can be estimated following the steps above,
Mxxxx=xxxxxxxx
P ((x - .)xx X/
(52)
After replacing the population conditional mean embedding as their empirical estimates, we can anrive at the solution.
19

---

Proposition 9 (GP-SHAP), Let the matrix A be defined as in Theorem 7 for the multivariate stochastic Shapley values can be estimated as,
The mean and covariance
Φ ( )) = (A [d]) (BXX),A(xx,,
(exB(x,|d)A)
(8)
where Bxx, [d) = [b .d..)...... b(x, [d]
Proof. The result follows directly from the previous proposition. Recall ((, :)) A for V the vector of stochastic payoffss for each coalition. To estimate the mean, we
Blo(0,((x ))] [
(53)
my (1
A
(54)
m, S2d)
b(x, S1) T ((()
A
(55)
b(x, S2u) (X)]
AB [d] mX)).
(56)
Recall V[vx]1, ((x, S:), (x, S,) b(), Kxxb(x, S,). the derivation for the covariance matrix then follows analogously as the derivation for the mean,
Vlo(0,(x, )) = AVMxAA
(57)
2",2a
= A b(x, 8:) Kxxb(x, 5T)
AT
(58)
i=15=1
= AB(x, [d] X(x Id)A
(59)
Proposition 10 (BayesSHAP and o2 follows:
20
Given the data generation above, the posterior distribution on Ø
Ze, f. D N(AVx (Z) WeZz))))''3
(11)
o² Zi,,,x, D Scaled-Inv- x lo + loo + (vx()x)
(12)
+
where l the number coalitions S {S we sample uniformly from 2 Ze is the binary matrix representing S, and W is the corresponding weight matrix, and A₁ = (Z, W. Ze) -' W the WLS matrix, Vx= = [ Si). (x St)] is the vector of deterministic payoffs, and
(W ( [Vv Z(AoZ T
(A.Vs) (ARA)
(13)
measures the average weighted error in the regression and the norm of the mean explanations.
Proof. See Slack et al.
20 Section. 3.1].
Proposition 11 (BayesGP-SHAP), Continuing from Propositions and 10 the posterior distribution of the stochastic Shapley values can be estimated using the Bayesian WLS approach as,
Φ 02 Zz.x,D~ (AB((x,S)) (m))m((X))))))) A,B(x
() + (ZZ)))
where o² is sampled from o² Ze ~ Scaled-Inv- x² (eo + fno
lott
Proof.
We drop the bar notation of Ø to unify notations. Given the posterior GP f D ~ ().,,
p(o 0² Ze D) / p( oZ D)r(f D)df
Ze,
(60)
20

---

Using a standard Gaussian conjugacy procedure, we can derive the variance as the sum of variances from GP-SHAP and BayesSHAP. While is possible integrate p(o² Z, f.x, x, D) with respect to the posterior, this leads to a complex scaled mixture of normals that is difficult to model. Instead, we construct scaled inverse chi-square distribution with s² (T[[vx)) which represents the error of the weighted regression with respect to the mean payoffs E v( We sample o² from the following distribution:
o² Zi, x, D ~ Scale-Inv-y' x² (6+ loo + (II(V))
(61)
to +
B.4 Proofs for section 4 on predictive explanation and Shapley prior
Proposition 12 (The Shapley prior over Φ). The prior J ~ GP(O. k) and the corresponding stochastic
gume (x S) = [ Xs = xs] induce vector-valued GP prior over the explanation
functions & ~ (P(O,, x) where X x X + Rdxd is the matrix-vals valued covariance kernel
<(x, x') A(x) A(x's A(x) A
(14)
where (x) [E[ X) Xs, rs........ B[k( X) XSyI = IS xs I] and the T sign refers to taking inner products in the RKHS of k.
Proof. The proof is similar to how we proved previous propositions but applied to prior GP f ~ GP(O, k)) instead. If we set,
(x = = Xs = xss,
(62)
then vy is GP on the joint space of data and coalitions with mean O, and covariance function,
cov (uy (x 8))) E'k(X, Xs xs, X's x'sl
(63)
Xx Xx|
(64)
Since Avx for Vx the vector of stochastic payoff from the game induced by the GP prior, the mean stays O, and the covariance is,
24 2d
x(x,x') = A T MXXXs, MXIXs, =x's
A
(65)
i=1,==
= =AV(x) (x
(66)
= = -Axx))(xx)),
(67)
therefore we have a matrix-valued covariance kernel x to build a prior over the induced Shapley values.
Proposition 13 (Predictive explanations as multi-output GPs). Given Do {xxi, d.))
(X. xx) where Φi Rd are the Shupley values computed under predictive model f and 4x (o1..... Tn)T the predictive explanations for new data x' is distributed as,
o(x') Do ~ N (x)))))), s(x', x') x(x', X)),, (x,, X))
(15)
where (x))x' = bx(xx' , X)T xx), bx(x'(X) = (xxx + 03 (3)))~n(XX,)) [XX, Xxx is the gram matrix for kernel K of size nd x nd, n(x', X) = [s(x', x1), x(xx, xn)) is size d x nd and o% is the noise parameter for regression.
Proof. Follows from standard vector-valued Gaussian process regression results. See Alvarez et al. 51 for a detailed discussion on regression with matrix- valued kernels.
Proposition 14 (Posterior mean as Shapley values for payoff vector Vx). The posterior mean
m (x)) corresponds to Shapley values for the payoff vector Vx', i.e., mx(x') Avx, where Vx En_ (xx') T (x,) a4 and as E Bdd is the [i,... + (d 1)] subvector of (Xxx
1) vec( x...
21

---

Proof. There are two ways to see this. First is by brute force and rearranging the terms in the posterior mean expression. The other is leverage the vector-valued representer theorem [22 and write the
posterior mean as,
m = Σ A(x:)aa, og E R
(68)
AV(x') (Xx))
(69)
Σ V(x)AA
(70)
A after some linear algebra exercises, we can see that aj is the li i + 1)] sub-vector of (xxx + vecc
(71)
C
Implementation details and further illustrations.
All illustrations run locally MacbookPro 2021 with Apple M1 pro chip.
C.1 Ablation study on different notions of uncertainties captured
To demonstrate the difference between the uncertainties captured by GP-SHAP, BayesSHAP, and BayesGP-SHAae SHAP, we utilise the California housing dataset [41]. This dataset was derived from the 1990 U.S. census, each observation represent a census block group. A block group is the smallest geographical unit for which the U.S. Census Bureau publishes sample data (a block group typically has population 600 to 3,000 people). The dataset includes 20640 instances with numerical features measuring the following:
MedIne: Median income in block group
HouseAge: Median house age in block group
AveRooms: Average number of rooms per household
AveBedrms: Average number of bedrooms per household
Population: Block group population
AveOccup: Average number houehold members
Latitude: Block group latitude
Longitude: Block group longitude
The target variable is the median house value for California districts, expressed in hundreds of thousands of dollars. In the following, we train a GP model and extract explanations using GP-SHAP, BayesSHAP, and BayesGP-SH.Ne for 4 different configurations:
1. trained on 25% of data, estimate the Shapley values using 50% of coalitions.
2. trained on 25% of data, estimate the Shapley values using 100% of coalitions.
3. trained on 100% data, estimate the Shapley values using 50% of coalitions.
4. trained 100% of data, estimate the Shapley values using 100% of coalitions.
To fit the GP model. we employ a sparse Variational GP approach with 200 learnable inducing point locations. The evidence lower bound is optimized using batch gradient descent with batch size of 64, a learning rate of 0.01, and 100 iterations. The RBF kernel with learnable bandwidths initialized using the median heuristic approach is used for the sparse GP. The inducing locations initialized using a standard clustering approach to obtain a representative set of inducing points.
After training the model, we reuse the learned inducing points and kernel bandwidths for the explanation algorithms. The explanations are obtained using the procedure described in Algorithm of our work.
22

---

92h
H4g
Modine
Population
Lattude
Avcoccup
Avwtrrup
Medine
Lingtude
Latiuus
Rolture
Poputation
Longitudo
Houswage
Hauseage
ArcRooms
Ar:Booms
AveBodrms
AvcRooms
0.00
0.25
0.50
0.73
125
1.50
1775
Moan stochastc Snapley
Standard deviatioo stochastic Shopley
(a) Mean each stochastic explanations
(b) Standard deviation each explanations
Figure 4: We plot the beeswarm plot of the mean and standard deviations of each stochastic explana- tions from BayesGP-SHHa fitted on the housing dataset. The features are ranked according to the distance span by the largest and smallest mean (std) stochastic Shapley values.
In Figure I of our paper, we present the stochastic Shapley values for the 1 11th observation, computed using the three explanation algorithms. The plot includes the 95% credible interval to visualize the uncertainties associated with the explanations.
Further illustration: In Figure 4 we plot the becswarm plot on the mean and standard deviation of each stochastic explanations respectively. We color the point based the relative size of the feature value compared to the rest. We see that in Figure 4a which plotted the mean stochastic shapley values for each observation, the relationship between most features' explanation the target variable is quite linear. For example, the higher the median income (MedInc), the more positive those feature contribute to predicting the respective median house value. On the other hand, Figure 44b illustrated the standard deviation of each stochastic explanations. In general, we see that the larger the feature values are, the more uncertain the explanation becomes. Nonetheless, feature contributing the most, defined as the feature having largest distance spanned by their most positive and most negative mean stochastic Shapley values, does not necessarily have the largest. variation respectively.
that the
C.2 Exploratory analysis of the stochastic explanations
For this illustration, we utilise the breast cancer dataset [4), containing 569 patients with 30 numeric features. They are computed from digitized image of fine needle aspirate (FNA) of a breast mass and describe characteristics of the cell nuclei present in the image:
radius (mean of distances from center points on the perimeter)
texture (standard deviation of gray-scale values)
perimeter
smoothness (local variation in radius lengths)
compactness (perimeter"
area-1
concavity (severity of concave portions of the contour)
concave points (number of concave portions of the contour)
symmetry
fractal dimension ("coastline approximation"- 1)
The goal is to predict whether tumour is malignant or benign. We first fit a GP model with RBF kernel using again the sparse Variational GP formulation with 200 learnable inducing locations. We initialise the inducing points using standard clustering techniques on the data. The evidence lower bound objective optimised with learning rate of 1e-4 and 1000 iterations using batch gradient descent of batch size 64. To obtain the explanations, run the BayesGP SHAP algorithm with 2
23

---

11g1
Hort TEAVE
werst fracial dimersion
frxctsl
symmetyy
paramater
Word OOTICEN 20:11
Socuue CIFOr
MOrST toxture
compactness eITor
ConiaNN points error
tractel dimens
concauity
symmetry
moan perfmoter
hrstt dimendoo enor
wenst compactness
Pastuan
Felture
mean concaity
concauity erter
symmesry
mean symmetry
mosn Focctll
worst radiax
mean dimersion
meen per permeter
pertmeter entor
maan compactness
workE symmetry
moan texture
-10
0,5
0.1
0.5
Shapley
Standard dectatiin
Shaplly
(a) Mean each stochastic explanations
(b) Standard deviation of each explanations
Figure 5: We plot the violin plot of the mean and standard deviations of each stochastic explanations from BayesGP- SHAP fitted on the breast cancer. The features are ranked according to the distance span by the largest and smallest mean (std) stochastic Shapley values.
number of coalitions. We do not compare GP-SHAP and BayesSHAP here because the BayesSHAP uncertainties have shrunk to almost zero, i.e., the mean standard deviations from the BayesSHAP uncertainties across all features and data is 0.0002. This reconfirms the fact from Slack et al. [20 that as we increase the sample size the estimation error goes to zero, thus the uncertainties from BayesSHAP goes to zero well. On the other hand, GP-SHAP uncertainties still remain valid because represents the GP predictive uncertainties, which do not shrink zero as we increase the number of coalitions we use to esitmate the SVs.
Further illustrations: In Figure 5 we plot two violin plots to illustrate the relationship between mean and standard deviation of the stochastic values with respect to the size of the original feature. see that the feature "worst fractal dimension" are the second most influential feature terms of mean stochastic explanations and also the feature that has highest uncertainty around its explanations. In comparison with the housing prediction problem illustrated in Figure the higher the feature value doesn't it necessary give higher uncertainty around its explanation.
C.3 Predictive explanations
For this illustration, we utilise the Diabetes dataset with 442 patient data and 10 numeric features measuring the following:
age: age in years
sex
bmi: body mass index
bp: average blood presuure
sl: total serum cholesterol
s2: low-density lipoproteins
s3: high-density lipoproteins
s4: total cholesterol
s5: Log of serum triglycerides level
24

---

* s6: blood sugar level
The experiment is to assess the effectiveness of the Shapley prior we proposed in predicting explana- tions estimated using SHAP algorithms for general models, including GP-SHAP, TreeSHAP, and DeepSHAP. We use the implementation of TreeSHAP and DeepSHAP from the shap package (2).
While algorithms such FastSHAP [22] also learn a vector-valued function that returns explanations given instances, the algorithm require access to the underlying model f during training while ours required previously computed explanations. Due to this importance difference in the problem setup, we do not compare the two algorithm.
We first generate three sets of explanations set up three regression problems:
1. Fit Gaussian process model and then GP-SHAP to obtain explanations.
2. Fit a random forest model and then run TreeSHAP to obtain explanations.
3. Fit a neural network model and then run DeepSHAP to obtain explanations.
After obtaining explanations as groundtruths for this experiment, we randomly divide 70% of them training data and 30% of them as testing data. We then do the following.
1. We fit a multi-output GP using the proposed Shapley prior on the training data and predict the explanations of the unseen test data.
2. We fit a multi-output random forest model on the training data and predict the explanations of the unseen test data.
3. We fit a multi-output neural network model on the training data and predict the explanations of the unseen test data.
We repeat this experiment 10 times using different seeds and compute the RMSE between the predicted and groundtruths explanations. The results are then plotted in Figure 3
C.4 Further ablation study: Impact of increased posterior prediction uncertainty on explanation uncertainties
In this ablation study, aim to examine the effect of increasing the uncertainty in posterior predictions on the corresponding uncertainty in stochastic Shapley values. To demonstrate this, we utilize the diabetic dataset [47) and split the data based on recorded sex. We train our GP model on the male data and employ BayesGP-SHAP to explain the prediction results for both the male training data and the female testing data. We adopt this split because we expect the biological characteristics between males and females to be distinct enough to treat the female data as out-of-sample data, thereby naturally resulting in increased predictive uncertainty for the female data. To further amplify this uncertainty, we multiply each instance in the female testing data by distortion factors of two and three, respectively, and assess the corresponding uncertainties in the explanations.
We begin by illustrating the relationship between the satpol-sampleness of the data and the corre- sponding increase in predictive posterior uncertainties. This is depicted in Figure 6aa where observe that as the data becomes more out-of-sample, the predictive uncertainties consistently rise. Even at distortion level which represents the original female data, can already observe increased uncertainties compared to the uncertainties derived from male data prediction.
Furthermore, these increased uncertainties in the predictive posterior are reflected in the associated feature explanations. This evident in Figure 6b) where we visualize the uncertainties associated with the feature explanations. For instance, the green bars representing the average uncertainties in explaining female data with no distortion are consistently larger than the red bars, which represent the average uncertainties of male data explanations. This observation aligns with the higher predictive uncertainties observed in Figure 6a for the female data compared to the male training data.
It is worth noting that the uncertainty for the feature "sex" remains consistently close zero. This is because the feature "sex" is constant within both the female and male datasets. As a result, it acts as the null player in each dataset and obtains almost Dirac zero as its stochastic Shapley value.
25

---

Female, distertion level =3.0
Female data, distortion level 3.0
distortion
Female data, distortion level 2.0
Fernae,, distortion level =1.0
of explonations
Female data, distortion level 1.0
Male (train)
Male data (train)
0.04
20-
0.2
0.00
45
standard deviations predictions
3100
5ml
(a) Predictive posterior standard deviation
(b) Mcan of standard deviations of explanations
Figure 6: Ablation study: (left) We begin by training Gaussian Process (GP) model the male data. We then make predictions using this trained model on both the male data and out-of-sample female data. To assess the impact of increasing posterior uncertainties, we multiply the female data by distortion levels of 1.0, 2.0, and 3.0. We visualize the results by plotting the density plot of the standard deviations obtained from the predictive posterior distributions. (right) Next, focus analyzing the average standard deviations of explanations per feature from the male and female data, considering different distortion levels. We observe that as we progressively increase the posterior uncertainties in the sample, these uncertainties are reflected in the uncertainties of the explanations provided.
Comparing FastSHAP with Shapley kernel approach over predicting test explanation
Average Absolute Shapley Values over test data
groundtruth explanations
0.30
shapley kernel approach (RMSE: 0.05)
fastshap (RMSE: 0.09)
0.25
0.15
0.10
age
5ex
bp
Figure 7: Comparing the predicted Shapley values with groundtruth from GP regression with Shapley kernel versus fast the diabetic dataset.
26

---

Distribution of kendall tau score compared with groundtruth explanation
Shapley kernel approach
FastSHAP
Density
3
-0.2
0.0
0.2
Kendall's Tau score
Figure 8: Comparing the Kendall's Tau score between the feature importance ranking from GP regression with Shapley kernel and FastSHAP with groundtruth feature importance ranking on the diabetic dataset.
D Comparison with FastSHAP 139
ran a predictive explanation experiment (see Figure and 7 using diabetic dataset from UCL, comparing both our method and FastSHAP. In particular, we use 70% of the data as training, i.e. we first compute the Shapley values for training data and treat them as labels, and predict the explanation for the remaining data. We evaluated performance using both RMSE and Kendall's 's tau distance (the closer to I the better). Notably, Kendall's tau distance assesses the consistency in feature importance ranking between predicted explanations and the true explanations. In the experiments, got RMSE 0.05, while FastSHAP got 0.09. In terms of Kendall's 's tau, our method got an average around 0.7 while FastSHAP got around 0.4, meaning the feature importance ranking recovered from our method is more aligned with the groundtruth SVs than FastSHAP. Our findings indicate that our method outperformed FastSHAP in this experiment.
method
Note however that interpretation of these comparisons should nuanced since the two techniques have different usecases: FastSHAP uses knowledge of the model, not the previously computed SVs, whereas GPSHAP uses previously computed SVs and agnostic to the type of model they came from.
27


Reviews:
[{"review_id": "xepCSuowVO", "reviewer": "Reviewer_1zjr", "summary": "The paper introduces a novel approach for explaining Gaussian processes (GPs) by utilizing the analytical covariance structure present in GPs based on extending the concept of Shapley values to stochastic cooperative games, resulting in explanations that are random variables. The GP explanations generated using this approach satisfy favorable axioms similar to standard Shapley values and possess a tractable covariance function across features and data observations. The proposed approach, GP-SHAP, and its variant BayesGP-SHAP, are shown to be effective in providing explanations for Gaussian process models.", "questions": "The font in Figure 2 is too small to read. ", "limitations": "Shapley-value based explanations have some limitations -- two references are provided as follows. \n* Shapley Residuals: Quantifying the limits of the Shapley value for explanations (https://proceedings.neurips.cc/paper/2021/file/dfc6aa246e88ab3e32caeaaecf433550-Paper.pdf)\n* Problems with Shapley-value-based explanations as feature importance measures (https://arxiv.org/pdf/2002.11097.pdf)\n\nThese limits are inherited in this proposed method. ", "rating": 7, "confidence": 3, "soundness": 4, "presentation": 4, "contribution": 4, "strengths": "The paper presented a novel and principled approach for explaining GP model using stochastic Shapley values. The proposed approach is sound and thoroughly evaluated. It also fills the void that there is no model-specific explanation method for GP. ", "weaknesses": "The paper should also discuss the limit of using Shapley value in explanation in general.  ", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "We thank the reviewer for sharing literature pointing out the limitations of Shapley-value based explanations. We have included a discussion and included the references in the discussion section in the camera ready version. We will also enlarge figure 2 in the camera ready version."}, {"review_id": "8j8aB1GpMj", "reviewer": "Reviewer_FXFV", "summary": "This works introduces the stochastic Shapley values to the explainable AI community. It formulates the stochastic Shapley values for GP models with a closed-form expression for the covariance matrix. Then, it proposes BayesGP-SHAP, which incorporates both sources (of uncertainty from the GP posterior and of the Shapley value estimation). Lastly, it constructs an inference model by proposing a GP prior to the Shapley value to solve the problem of predictive explanation.\n", "questions": "I do not have any questions for the authors.\n", "limitations": "The authors have adequately stated the limitations in the appendix.", "rating": 8, "confidence": 3, "soundness": 4, "presentation": 3, "contribution": 4, "strengths": "I enjoy reading the paper as it contains several major strengths: \n+ The motivation of predictive uncertainty in trustworthy ML models is convincing, which paves the way to the proposed GP-SHAP and BayesGP-SHAP.\n+ The literature review is sufficient to highlight the novelty of the approach.\n+ It introduces concepts of stochastic cooperative games and stochastic Shapley values to the explainable AI community.\n+ GP explanations not only satisfy favorable axioms to standard Shapley values but are also able to measure explanation uncertainties and determine statistical dependencies between explanations. \n+ The proposed BayesGP-SHAP integrates both sources of uncertainty from the GP posterior and the estimation of the Shapley value.\n+ It proposes the Shapley prior to dealing with the problem of predictive explanation, which is a nice Bayesian treatment to the predictive explanation problem.\n+ In the experiments, the paper introduces several exploratory analysis methods for practitioners to understand stochastic explanations.\n", "weaknesses": "The paper may elaborate on other use cases of predictive explanations. Anyways, I find the Bayesian optimization application in the conclusion interesting.\n\nThere are very minor issues in the references as some words are not capitalized, e.g., Shapley.\n", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "We thank the reviewer for their positive feedback and pointing out the typos in the references, we will amend them accordingly in the camera ready version. "}, {"review_id": "WjvvJ03H7o", "reviewer": "Reviewer_hgou", "summary": "The authors consider Shapley values computed from games $\\nu$ modelled with a Gaussian process. Since the game is stochastic, the Shapley values are now stochastic, and the stochasticity reflects epistemic uncertainty in the explanations. This is unlike previous methods like BayesSHAP, where the uncertainty is estimation error. The approach of the authors, called GP-SHAP, can be used to explain predictions of Gaussian processes. It can also be used to construct a prior over Shapley values, which in turn can be used to perform regression on Shapley values. The authors demonstrate the utility of their approach is three experiments.", "questions": "## Conclusion\n\nI think this is a solid submission without any major shortcomings, which is why I am giving an accept.\n\nIf at all possible, I would like to see the following edits in a revision of the authors:\n\n* If I'm right that Theorem 4 follows immediately from the deterministic case, I would really just do that and omit the current proof.\n\n* Please start out Section 3 by explaning what you are working towards: pushing a GP through the value allocation computation\n\n* Please clarify the statement of Proposition 11 (see above).\n\n* Please clarify the transpose in Proposition 12 (see above).\n\n* In Figure 2.(b), please analyse the difference in ordering between the mean of absolute SSVs and the absolute values of SSVs and conclude which of the two orderings is more sensible.\n\n* In Figure 2.(d), please argue whether the graphical model is sensible or not.\n\nEDIT\n\nI have read the authors' rebuttal and the other reviews and remain at my current assessment.", "limitations": "See above.", "rating": 6, "confidence": 3, "soundness": 4, "presentation": 3, "contribution": 3, "strengths": "I would like to thank the authors for their submission. I think this is interesting work.\n\n## Strengths\n\n* The paper is well written and mostly clear. The presentation is reasonably good.\n\n* Pushing a GP through the value allocation computation to explain stochastic predictions seems like a very sensible thing to do.\n\n* Doing the above to define a prior over value allocations and then using this prior for regression of Shapley values is, to the best of my knowledge, novel. I have not seen this before. I should, however, say that I'm not at all familiar with this part of the literature.\n\n* The stochasticity in explanations by GP-SHAP represents epistemic uncertainty in the explanation, whereas the stochasticity in explanations by BayesSHAP is due to estimation error. I think this is an interesting finding.\n\n* I have either skimmed through or looked in more detail at the proofs of the mathematical statements. Apart from some clarity issues (see below), I think everything generally checks out.\n\n* The experiments nicely illustrate the benefits of GP-SHAP.\n\n\n## General Remarks\n\n* In Proposition 11, I was initially very confused by that the expectation of $v_x$ appears. After reading the proof, it became clear that this is an approximation that simplifies the computation, which should actually integrate over $f$. I think this can be better explained in the statement of Propostion 11.\n\n* In Proposition 12, if I'm not mistaken, the transpose in the expression for $\\kappa$ refers to taking inner products in the RKHS of $k$. This is _definitely_ not clear at all and should be clarified!\n\n* In Figure 1, how does BayesSHAP produce its explanations? It is applied to the predictive mean of the GP prediction?", "weaknesses": "## Weaknesses\n\n* I think that providing a proof for Theorem 4 is unnecessarily complicating the exposition, as the Shapley's original proof for the deterministic case can be applied to immediately prove Theorem 4: Let $(\\Omega, \\mathcal{F}, \\mathbb{P})$ be the underlying probability space. For any fixed sample $\\omega \\in \\Omega$, $\\nu$ is d-game (using the authors' language). Therefore, according to Shapley's original proof, for that $\\omega \\in \\Omega$, it can be written in the form of (1). Since (1) holds for all $\\omega \\in \\Omega$, it obviously holds as an equality in terms of random variables.\n\n* Generally, I think that the exposition takes a long time to arrive at the key idea of the paper: model $\\nu$ with a GP, and push it through the value allocation computation in (1). I think that the exposition would be improved by having Section 3 starting out with immediately explaning this key idea. A paper is not supposed to be a novel: please signpost important ideas and conclusions as much as possible.\n\n* In lines 326-339, you highlight the difference between the mean of absolute SSVs and the absolute values of mean SSVs, and state that this gives a different ordering for which feature is most influential. However, you do not analyse whether this different ordering is better or worse, which means that it is not clear whether, in this case, GP-SHAP produced a better ordering or not.\n\n* In lines 340-346, you produce a local explanation graphical model. However, you do not analyse whether this graphical model is reasonable, which means that the reader is not sure whether GP-SHAP did something sensible or not.\n\n* Building a GP prior over Shapley values to do regression is an interesting idea. However, the presentation would be more convincing if you were to point out some existing applications that would actually want to do regression of Shapley values.\n", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "We thank the reviewer for their in-depth suggestions. We would like to address their concerns below:\n\n### Proof of Theorem 4:\nWe agree with the reviewer that using the d-game to establish stochastic Shapley values as a push forward measure is easier to understand. In fact, this approach was the one adopted by Ma et al. 2008, where they first proposed a potential solution, then verified its alignment with the three (stochastic) axioms and its uniqueness—a top-down proof strategy. For our study, as explained in lines 561-565, we aimed to offer a contrasting perspective that hasn't appeared in the literature, mirroring Shapley's original bottom-up derivation. This means we began with the axioms and subsequently determined the unique solution. Recognising the potential for confusion, we have incorporated the following explanation in the main text under theorem 4 to clarify our approach:\n- **``We emphasise that this result has been proven in Ma et al. using a top-down approach, i.e. starting with Equation (1) and then verifying that it satisfies the stochastic axioms and uniqueness. In the appendix, we offer a contrasting perspective where we mirror Shapley's original bottom-up derivation, i.e. we began with the stochastic axioms and subsequently determined the unique solution''**\n\n[Ma et al. 2008] Ying Ma, Zuofeng Gao, Wei Li, Ning Jiang, Lei Guo, et al. The shapley value for stochastic\ncooperative game. Modern Applied Science, 2(4):1–76, 2008.\n\n---\n\n### Clarity\nWe have incorporated the reviewer’s suggestion and included a high level description of key ideas in the beginning of section 3. Both propositions 11 and 12 are now clarified as well. Particularly, we wrote:\n- High level description: **\"This section is dedicated to introducing GP-SHAP, a method which efficiently computes stochastic Shapley values for stochastic games where the value function is computed using conditional expectations of Gaussian processes.\"**\n- Proposition 11: **\"We note that in the above proposition, instead of integrating $p(\\sigma^2\\mid \\mathbf{Z}_\\ell, f, \\mathbf{x}, \\mathbf{D})$ with respect to the posterior GP, which leads to a complex scaled mixture of Gaussians, we simplify the computation and instead construct a scaled inverse chi-square distribution with $s^2(\\mathbb{E}[\\mathbf{v}_x])$ instead, which represents the error of the weighted regression with respect to the mean payoffs.\"**\n- Proposition 12: We will add this sentence to the end of prop 12, **\".. and the inner product is taken over the RKHS of $k$\"**\n\n---\n\n### Global explanation\nWe believe that it's a more principled approach to measure global importance by computing the mean of absolute SSVs rather than taking the absolute values of the mean SSVs, which neglects uncertainty in SSVs. Practitioners employing GP-based models would likely prioritise leveraging quantified uncertainty in their explanations. By simply using the absolute values of mean SSVs, they neglect to incorporate this valuable uncertainty information. We have added the following discussion in the main paper to make the exposition clearer:\n- **\"We believe the former (mean of absolute SSV) approach is more appropriate than the latter (absolute of mean SSVs) as practitioners employing Gaussian processes models would want their global explanations to take the quantified uncertainty into account.\"**\n\n---\n### Graphical model:\nWe emphasise that, analogously to the Shapley values as explanation tool, the graphical model depicting dependence structure of stochastic Shapley values is a tool for explaining the fitted GP model, not the true underlying data generating processes, so whether or not this graphical model is sensible will depend on whether the fitted GP model is itself sensible. Nonetheless, Gaussian graphical model is a faithful representation of the multivariate normal distribution of the Shapley values that is induced by the fitted GP model, and such representation of multivariate normal distributions is standard in the literature.\n\n---\n\nWe thank the reviewer again for their suggestions on the revision. As we have answered and clarified the reviewer's comments, we hope that the reviewer could revise their score from a weak accept to accept, in light of our responses and other reviewers' positive feedback. Thank you."}, {"review_id": "EHpD9qtXnc", "reviewer": "Reviewer_uH3p", "summary": "The work proposes a new kind of Shapley value (SV) based explanation method. The new methods GP-SHAP and BayesGP-SHAP are like TreeSHAP, LinearSHAP or RKHS-SHAP model-specific SV methods for Gaussian Process (GP) models. The authors further propose another approach of utilizing GP models to amortize (i.e learn) SVs for a model and data point in a predictive task (like FastSHAP). The authors name this the Shapley prior. The two main contributions (model-specific) SV method and the Shapley prior is illustrated in experiments with well-established benchmark datasets.", "questions": "1. How does your approach compare to FastSHAP in terms of predictive performance? What's the tradeoff in terms of predictive performance and explainability?\n\n2. Can you say anything about what parts of uncertainty are captured in the explanations uncertainty (i.e. aleatoric or epistemic)?", "limitations": "The work includes a limitation section discussing (a) limitations and (b) future work. However, the section is quite short and focuses solely on future work and avenues to explore (e.g. extension to Bayesian NN). A more detailed discussion of the limitations of this particular result would be beneficial.", "rating": 7, "confidence": 4, "soundness": 4, "presentation": 4, "contribution": 3, "strengths": "The paper has many strengths. The work is theoretically well-formulated and motivated. It combines uncertainty research with XAI research to enrich the quality of explanations. This is very important research. The contribution is well formulated and well motivated (formally and verbally). The experimental analysis of the model-specific GP-SHAP and BayesGP-SHAP is convincing. While the work is very theoretical in nature, it is still very pleasant to read, which is hard to achieve, since bridging the gap between two streams of research is not easy. \nThe work adheres perfectly to the style guidelines, includes a well-rounded appendix with code and a small readme file for easier reproducibility. The proofs are to the most part understandable.", "weaknesses": "The main weaknesses I identified are elaborated below:\n\n1. The predictive (FastSHAP-like) approach to learning SVs with GP models and the Shapley prior feels a bit tagged on and underdeveloped. In lines 256-258 the argument is made that GP models are trained only with x (instance) and y (SV for instance) pairs and the model is not required to train. This is posed as a positive side of this Shapley prior approach. This, of course, is correct that you only need the SV to train, however, doesn't this make the GP models less effective in comparison to, for instance, FastSHAP which is trained on this Shapley loss function (i.e. MSE over model output of coalition with feature i than without feature i). The main benefit I see is that you now have a GP model rather than a big black-box NN like FastSHAP proposes. The GP models, therein, are more explainable which is a desideratum in my opinion for an explanation system (this \"explaining a model with another model that needs to be explained in itself\" can be spun ad absurdum). However, the paper is missing such an experimental comparison of the Shapley prior with FastSHAP or an in-depth discussion of the benefits in my opinion. I want to reiterate that this is not a big problem, however, it is left hanging.\n\n2. The work has some accessibility issues from a ML perspective. Some concepts are hard to follow even with a good background in SV and statistics. I recon that ML readers might have a hard time understanding the contribution well because they already have problems with the problem definitions (beginning of 2). If space allows I'd suggest to certainly enrich the definition section (lines 81-92) with some explanations.\n\n3. The connection to uncertainty research could be a bit stronger. I.e. it would be nice to include a more extensive discussion of what parts of uncertainty are captured in the confidence intervals (i.e. aleatoric or epistemic)?\n\nMinor things:\n- Line 11 ... \"extensive illustrations\" sounds a bit counterintuitive consider changing this\n- Write out the main contribution again in the introduction in a seperate paragraph that is easily digestible for ML engineers\n- Line 134 \"rvs\" is introduced but (I think) not used. Even if it is used consider removing this abbreviation\n- Lines 235-240: I would like to see a more extensive discussion here.", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "We thank the reviewer for their constructive feedbacks and we will answer their questions below:\n\n## 1. How does your approach compare to FastSHAP in terms of predictive performance? What's the tradeoff in terms of predictive performance and explainability?\n- To answer this, we have ran a predictive explanation experiment (please see attached pdf in global response) using the diabetic dataset from UCI, comparing both our method and FastSHAP. The results are included in the one-page pdf submission. In particular, we use 70% of the data as training, i.e. we first compute the Shapley values for training data and treat them as labels, and predict the explanation for the remaining data. We evaluated performance using both RMSE and Kendall’s tau distance (the closer to 1 the better). Notably, Kendall’s tau distance assesses the consistency in feature importance ranking between predicted explanations and the true explanations. In the experiments, our method got RMSE 0.05, while FastSHAP got 0.09. In terms of Kendall's tau, our method got an average around 0.7 while FastSHAP got around 0.4, meaning the feature importance ranking recovered from our method is more aligned with the groundtruth SVs than FastSHAP. Our findings indicate that our method outperformed FastSHAP in this experiment.\n- **Note however that interpretation of these comparisons should be nuanced since the two techniques have different usecases: FastSHAP uses knowledge of the model, not the previously computed SVs, whereas GPSHAP uses previously computed SVs and is agnostic to the type of model they came from.**\n---\n## 2. Can you say anything about what parts of uncertainty are captured in the explanations’ uncertainty (i.e. aleatoric or epistemic)?\n- Both GPSHAP and BayesSHAP are capturing epistemic uncertainties, but from different sources. GPSHAP captures epistemic uncertainty from fitting the GP function f, while BayesSHAP captures epistemic uncertainty arising when estimating the Shapley values using a Bayesian weighted least square formulation (for a fixed function f). BayesGPSHAP then captures both uncertainties.\n\n"}]

Author rebuttal:
We include further experiments to compare FastSHAP with GP regression using the Shapley prior in light of reviewer uH3p's comments. 

Meta-review: There is a clear consensus among the reviewers that this paper makes an interesting and substantial contribution. From the beginning, they have all been on the positive side, although two reviewers also raised a couple of concerns and made suggestions for improvement. Nevertheless, after the rebuttal, all reviewers are still in favor of acceptance, and so am I as well.
Decision: Accept (spotlight)