--- title: 'Permutation SHAP: Sampling-Based Shapley Estimation' url: https://www.emergentmind.com/topics/permutation-shap type: topic --- # Permutation SHAP: Sampling-Based Shapley Estimation Searching arXiv for recent and foundational papers relevant to Permutation SHAP. Permutation SHAP denotes a family of Shapley-value procedures in which permutations play the central computational role. In its standard usage, PermutationSHAP is a simple, model-agnostic way to compute Shapley attributions by averaging marginal contributions over random feature orderings; in its paired-sampling form, each sampled ordering is coupled with its reverse, yielding asymptotic normality, exactness for purely second-order interactions, and an additive-recovery guarantee for block-additive value functions [2508.12947]. In a separate recent usage, “permutation SHAP” refers to computing SHAP or KernelSHAP on data obtained by independently permuting each column, thereby sampling from the product of marginals $\mu^*$ and enabling theoretically justified feature discarding from aggregate scores [2503.23111]. ## 1. Definition, notation, and scope Let $f:\mathbb{R}^p\to\mathbb{R}$ be a predictor and let $x\in\mathbb{R}^p$ be the instance to explain. For any subset $S\subseteq\{1,\dots,p\}$, write $x_S$ for the input where only features in $S$ are retained and the rest are masked, intervened, or set to some baseline. If $\pi$ is a permutation of $\{1,\dots,p\}$ and $S_i(\pi)$ is the set of features preceding $i$ in $\pi$, then the classical Shapley value is $$ \phi_i(f,x)=\frac{1}{p!}\sum_{\pi\in S_p}\Big[f(x_{S_i(\pi)\cup\{i\}})-f(x_{S_i(\pi)})\Big]. $$ Equivalently, with $v(S)=f(x_S)$, $$ \phi_i=\frac{1}{p!}\sum_{\pi}\Big[v(S_i(\pi)\cup\{i\})-v(S_i(\pi))\Big]. $$ Because $p!$ becomes infeasible once $p\gtrsim 10$, PermutationSHAP replaces the full average over $S_p$ by sampling random permutations [2508.12947]. The permutation formulation is also written as an expectation over $S_d$ in the sampling literature. With $S_\pi(i)=\{j:\pi^{-1}(j)<\pi^{-1}(i)\}$, one has $$ \phi_i=\frac{1}{d!}\sum_{\pi\in S_d}\Big[f(S_\pi(i)\cup\{i\})-f(S_\pi(i))\Big], $$ which makes the approximation problem a numerical integration problem on permutation space. In that setting, a simple unbiased Monte Carlo estimator converges at $O(n^{-1/2})$ and requires $O(n\cdot d)$ calls to $f$ if one exploits the sequential structure of $\pi$ [2104.12199]. The cited literature uses the phrase “permutation SHAP” in two related but distinct senses: | Usage | What is permuted | Primary purpose | |---|---|---| | PermutationSHAP | Feature orderings $\pi\in S_p$ | Local Shapley estimation by averaging marginal gains | | Extended-support permutation SHAP | Columns of the data matrix | Global aggregation over $\mu^*$ for safe feature discarding | The first sense concerns local attribution for a fixed instance $x$. The second concerns aggregate scores over an extended distribution and is therefore conceptually closer to global feature selection than to pointwise explanation. ## 2. Paired-sampling PermutationSHAP estimator A paired-sampling trick reduces variance by coupling each sampled permutation $\pi$ with its reverse $\rho(\pi)$. If $\pi^1,\dots,\pi^m\sim \mathrm{Uniform}(S_p)$ are i.i.d., define $$ \Delta^{(k)}(i)=v(S_i(\pi^k)\cup\{i\})-v(S_i(\pi^k)), $$ and $$ \Delta'^{(k)}(i)=v(S_i(\rho(\pi^k))\cup\{i\})-v(S_i(\rho(\pi^k))). $$ The paired-sampling PermutationSHAP estimate is $$ \hat\phi_i^{\mathrm{perm}}=\frac{1}{2m}\sum_{k=1}^m\Big[\Delta^{(k)}(i)+\Delta'^{(k)}(i)\Big]. $$ Its implementation is direct: initialize $\hat\phi_i=0$ for all features, repeatedly draw a uniform permutation, compute the marginal gains for that ordering and for its reverse, and accumulate their average. Because each iteration uses exactly $2p$ model-evaluations, the total cost is approximately $2pm$ [2508.12947]. This estimator is model-agnostic in the sense that it only requires the ability to evaluate the value function $v(S)=f(x_S)$ for selected subsets. Its computational profile is therefore dominated by the number of model evaluations induced by the sampled orderings rather than by model-specific structure. In the summary provided for the method, the cost per sample is $O(p)$ model-evaluations [2508.12947]. The ordering-based construction is important. Unlike regression-based approximations, it accumulates incremental gains along a path through feature space, and the reverse permutation introduces a symmetric counterpart to that path. This structure is what underlies the variance reduction and the exact recovery results stated below. ## 3. Asymptotic normality, exactness, and additive recovery The paired estimator admits a central limit theorem. Define the random-vector contribution $B_\pi\in\mathbb{R}^p$ by $$ B_\pi(i)=v(S_i(\pi)\cup\{i\})-v(S_i(\pi)), $$ and let $$ \Sigma=\frac14\operatorname{Var}_{\pi}\big[B_\pi+B_{\rho(\pi)}\big]\in\mathbb{R}^{p\times p}. $$ Then, as $m\to\infty$, $$ \sqrt{m}\big(\hat\phi^{\mathrm{perm}}-\phi\big)\Rightarrow \mathcal{N}(0,\Sigma). $$ In particular, each component $\hat\phi_i$ converges at rate $O_p(1/\sqrt{m})$, and $\Sigma$ can be estimated empirically from the sampled permutations [2508.12947]. A stronger result holds for purely second-order interactions. If the value function is a bilinear form, $$ v(S)=x_S^\top A x_S, $$ so that only feature-pairs interact, then even $m=1$ suffices: for any single permutation $\pi$, $$ \frac{B_\pi(i)+B_{\rho(\pi)}(i)}{2}=\phi_i $$ exactly. The proof sketch given for the result shows that, in the purely quadratic case, the paired sum becomes independent of $\pi$ [2508.12947]. The distinctive structural guarantee is the additive-recovery property. Suppose $$ v(S)=\sum_{k=1}^K v_k(S\cap A_k), $$ where $A_1\cup\cdots\cup A_K=\{1,\dots,p\}$ is a partition into disjoint additive blocks. Then for any single permutation $\pi$, the paired-sampling estimator with one draw satisfies $$ \sum_{i\in A_k}\hat\phi_i^{(1)}=\sum_{i\in A_k}\phi_i,\qquad k=1,\dots,K. $$ The proof sketch proceeds by telescoping the incremental gains along $\pi$ within each block while features from other blocks act as dummy variables [2508.12947]. These theorems sharply distinguish errors internal to a block from errors across blocks. The block totals are exact after one paired sample, whereas finer-grain errors within a block shrink at rate $O_p(1/\sqrt{m})$. This suggests that paired PermutationSHAP is particularly well matched to models with low-order interactions or a block-additive structure. ## 4. Comparison with KernelSHAP KernelSHAP is presented in the cited material as a global weighted-least-squares procedure. In the paired-sampling comparison, it solves $$ \min_{\phi}\; \mathbb{E}_{Z\sim p_w}\Big[\big(v(Z)-\sum_i Z_i\phi_i\big)^2\Big] \quad\text{subject to}\quad \sum_i \phi_i=v(\{1,\dots,p\}), $$ and its sampling version is a regression of $v(Z)$ onto $Z$ [2508.12947]. The two estimators share some asymptotic features. KernelSHAP is unbiased and converges at $O(1/\sqrt{n})$, while paired PermutationSHAP satisfies a CLT at rate $O_p(1/\sqrt{m})$ with explicitly computable covariance. The difference emphasized in the comparison is not consistency but structure: PermutationSHAP’s incremental-gain sampling respects the ordering-based telescoping that recovers block sums exactly, whereas the regression-based KernelSHAP “mixes” contributions across blocks [2508.12947]. As a consequence, KernelSHAP does not, in general, return the correct aggregate $\sum_{i\in A_k}\hat\phi_i$ over an additive partition, even though it is unbiased in the large-sample limit. The additive-recovery guarantee is therefore exclusive to paired-sampling PermutationSHAP in the formulation summarized by the paper. In practice, the same summary states that paired-sampling halves variance, expressed as $\Sigma$ versus $\tfrac14\operatorname{Var}[B_\pi]$, and that for purely second-order interactions even one pair yields exact $\phi$ [2508.12947]. The comparison is methodologically important because it identifies a property that is not a generic consequence of Shapley-value axioms alone. It arises from the specific pathwise estimator defined over permutations. ## 5. Sampling design on permutation space Once Shapley estimation is written as an expectation over $S_d$, the quality of the permutation samples becomes a numerical integration problem. Standard Monte Carlo sampling is unbiased but can exhibit slow convergence, and variance can be high for models with strong feature interactions. Variance-reduction schemes such as antithetic, stratified, and Owen’s multilinear improve constants but not the $O(n^{-1/2})$ rate [2104.12199]. One line of work formulates the problem in an RKHS over permutations. With a positive-definite kernel $k:S_d\times S_d\to\mathbb{R}$, permutations are embedded into an RKHS $\mathcal{H}$. The cited paper lists three kernels: the Kendall kernel, the universal Mallows kernel $$ k_M^\lambda(\pi,\sigma)=\exp[-\lambda\cdot \tau_K(\pi,\sigma)], $$ and the Spearman inner-product kernel. In all three cases, $\mathbb{E}_{\sigma\sim \mathrm{Uniform}}\,k(\pi,\sigma)$ can be computed in closed form [2104.12199]. Kernel herding greedily builds an unweighted set of permutations to reduce RKHS discrepancy, and for universal kernels such as Mallows it achieves $O(1/n)$ convergence in RKHS norm. Sequential Bayesian quadrature models the integrand on permutations as a zero-mean Gaussian process and selects permutations by minimizing posterior variance of the integral estimate. The worst-case quadrature error satisfies $$ \left|\mathbb{E}_{\pi\sim U}[f(\pi)]-\sum_{\pi_t\in \Pi} w_t f(\pi_t)\right| \le D(\Pi,w)\,\|f\|_{\mathcal H}, $$ where $D(\Pi,w)$ is the kernel discrepancy; herding drives $D$ down for equal weights, and SBQ optimally chooses weights to minimize it [2104.12199]. A separate construction exploits the geometry of the permutohedron and the hypersphere $\mathbb{S}^{d-2}$. After centering and scaling, permutations lie on a translated copy of $\mathbb{S}^{d-2}$, and a point on the sphere can be mapped back to a permutation via rank ordering. The paper gives two fast schemes: orthogonal spherical codes, which produce $2(d-1)$ permutations with negative correlation and good coverage, and Sobol-on-the-sphere sampling, which uses a low-discrepancy sequence, inverse-CDF transforms for hyperspherical coordinates, and an $\mathrm{argsort}$ map back to permutation space. Both yield uniform $\pi\sim S_d$ in expectation and run in $O(n\cdot d)$ time [2104.12199]. Empirically, the same study reports that on small–medium $d$ with $d\lesssim 50$, kernel herding and SBQ with Mallows $\lambda\approx 4$ achieve $O(1/n)$ error decay and outperform MC’s $O(1/\sqrt{n})$, while on larger $d$ or very expensive models the linear-time sphere methods often match or beat antithetic MC by stronger negative correlation or low discrepancy. The experiments show convergence to a smaller RMSE in the same number of model evaluations [2104.12199]. ## 6. Extended-support permutation SHAP and feature discarding A separate line of work studies aggregate SHAP values for feature removal and shows that small aggregate SHAP values on the original data support do not necessarily imply that the corresponding feature does not affect the function. The issue is that SHAP evaluates $f$ on points outside the data support, where dependence on a feature can be masked. The proposed remedy is to aggregate SHAP values over the extended distribution $$ D_{\mathrm{ext}}\equiv \mu^*=\mu_1\times\mu_2\times\cdots\times \mu_d, $$ whose support is $\operatorname{supp}(\mu_1)\times\cdots\times \operatorname{supp}(\mu_d)$ [2503.23111]. Sampling from $\mu^*$ can be implemented by independently permuting each column of the observed data matrix. If $X\in\mathbb{R}^{n\times d}$ is the original data and $\sigma_j$ is an independent random permutation of $\{1,\dots,n\}$ for each column $j$, define $$ X^*_{kj}=X_{\sigma_j(k),j}. $$ Each row of $X^*$ is then an approximate sample from $\mu^*$ because coordinates are independent across columns. The paper defines the extended-support aggregate $$ \Phi_i^{\mathrm{ext}}(f)=\mathbb{E}_{x\sim\mu^*}\big[|\phi_i(\mu,f,x)|\big]. $$ Its main theorem states that $$ f \text{ depends only on } [d]\setminus\{i\} \text{ over } \operatorname{supp}(\mu^*) \iff \phi_i(\mu,f,x)=0 \quad \forall x\in \operatorname{supp}(\mu^*). $$ A robust version states that if $f\in[0,1]$ and $\Phi_i^{\mathrm{ext}}(f)\le \epsilon$, then there exists $g$ not depending on $x_i$ such that $$ \|f-g\|_{L^2(\mu^*)}^2