ControlSHAP: Variance Reduction for Shapley Values

• ControlSHAP is a methodology that stabilizes Monte Carlo Shapley value estimates for complex machine learning models using classical control variate techniques.
• It employs a Taylor-based surrogate model to dramatically reduce estimation variance while maintaining theoretical guarantees of unbiasedness.
• ControlSHAP integrates with existing workflows across various model classes and data distributions to enhance ranking consistency and efficiency in high-dimensional settings.

ControlSHAP is a methodology designed to stabilize Monte Carlo estimates of Shapley values for complex black-box machine learning models. Leveraging the classical control variates technique from Monte Carlo estimation, ControlSHAP exploits a computationally efficient Taylor-based surrogate model to dramatically reduce the variance of Shapley value approximations while maintaining theoretical guarantees of unbiasedness. Its design requires minimal modification to existing Shapley estimation workflows and is broadly applicable across model classes and datasets, offering a scalable solution to the inherent stochasticity of Shapley sampling in high-dimensional settings (Goldwasser et al., 2023).

1. Shapley Values and Monte Carlo Approximation

Given a predictive model $f:\mathbb{R}^d \to \mathbb{R}$ and a query point $x \in \mathbb{R}^d$, the Shapley value $\phi_j(x)$ for feature $j$ quantifies its average additive contribution to the model output, marginalizing over all possible feature orderings. The conditional expectation

$v_x(S) := \mathbb{E}[f(X)\mid X_S = x_S]$

captures the expected prediction when features in $S$ are known. The exact Shapley value is

$$\phi_j(x) = \frac{1}{d}\sum_{S \subseteq [d] \setminus \{j\}} \binom{d-1}{|S|}^{-1} [v_x(S \cup \{j\}) - v_x(S)]$$

or, equivalently, as an expectation over a uniformly random feature permutation. Due to $2^d$ combinatorics, practical estimation relies on Monte Carlo:

$$\widehat{\phi}_j^{\mathrm{MC}}(x) = \frac{1}{M}\sum_{m=1}^M [v_x(S^m \cup \{j\}) - v_x(S^m)]$$

where $S^m$ is formed by randomly permuting $x \in \mathbb{R}^d$0 and taking the prefix up to $x \in \mathbb{R}^d$1.

2. Sources of Variability in Shapley Estimates

Monte Carlo Shapley approximations $x \in \mathbb{R}^d$2 are subject to sampling error, as they average over only $x \in \mathbb{R}^d$3 random coalitions instead of the full $x \in \mathbb{R}^d$4 set. With $x \in \mathbb{R}^d$5,

$x \in \mathbb{R}^d$6

Further variance is induced when $x \in \mathbb{R}^d$7 is computed via nested Monte Carlo (i.e., expectation over unknown features is itself approximated).

3. Control Variates for Variance Reduction

Classical control variates use a secondary estimator $x \in \mathbb{R}^d$8 with known expectation $x \in \mathbb{R}^d$9 to reduce the variance of a primary Monte Carlo estimate. For estimators $\phi_j(x)$0 and $\phi_j(x)$1, one constructs:

$\phi_j(x)$2

where $\phi_j(x)$3 remains unbiased for $\phi_j(x)$4. The optimal $\phi_j(x)$5 minimizing variance is

$\phi_j(x)$6

giving residual variance $\phi_j(x)$7, where $\phi_j(x)$8 is the correlation between $\phi_j(x)$9 and $j$0.

4. ControlSHAP: Surrogate-based Control Variate Selection

ControlSHAP instantiates the control variate framework by setting

• $j$1: the empirical Shapley value for model $j$2,
• $j$3: the empirical Shapley value for a surrogate model $j$4 (a local Taylor expansion),
• $j$5: the analytically computed Shapley value for $j$6.

The surrogate $j$7 exploits the local differentiability of $j$8 and the data distribution moments.

4.1 Independent Features Case

Assuming feature independence, with data mean $j$9, covariance $v_x(S) := \mathbb{E}[f(X)\mid X_S = x_S]$0, Jacobian $v_x(S) := \mathbb{E}[f(X)\mid X_S = x_S]$1, and Hessian $v_x(S) := \mathbb{E}[f(X)\mid X_S = x_S]$2, the closed-form surrogate Shapley value (Goldwasser & Hooker Thm 1) is

$v_x(S) := \mathbb{E}[f(X)\mid X_S = x_S]$3

4.2 Correlated Features Case

For correlated features, assuming $v_x(S) := \mathbb{E}[f(X)\mid X_S = x_S]$4 is approximately Gaussian,

$v_x(S) := \mathbb{E}[f(X)\mid X_S = x_S]$5

where $v_x(S) := \mathbb{E}[f(X)\mid X_S = x_S]$6 is precomputed once (at $v_x(S) := \mathbb{E}[f(X)\mid X_S = x_S]$7 complexity) and then reused.

5. Algorithmic Workflow for ControlSHAP

The principal steps in ControlSHAP are as follows:

1. Precompute $v_x(S) := \mathbb{E}[f(X)\mid X_S = x_S]$8 matrices (correlated case).
2. For $v_x(S) := \mathbb{E}[f(X)\mid X_S = x_S]$9:
   • Sample a random coalition $S$0.
   • Evaluate $S$1 on $S$2 and $S$3 to get $S$4.
   • Evaluate the Taylor surrogate $S$5 on analogous inputs for $S$6.
3. Calculate Monte Carlo estimates for each $S$7 for both $S$8 and $S$9.
4. Compute closed-form $$\phi_j(x) = \frac{1}{d}\sum_{S \subseteq [d] \setminus \{j\}} \binom{d-1}{|S|}^{-1} [v_x(S \cup \{j\}) - v_x(S)]$$0.
5. Estimate variance and covariance terms empirically, by regression theory, or bootstrap.
6. Set $$\phi_j(x) = \frac{1}{d}\sum_{S \subseteq [d] \setminus \{j\}} \binom{d-1}{|S|}^{-1} [v_x(S \cup \{j\}) - v_x(S)]$$1.
7. Form the final ControlSHAP estimator:

$$\phi_j(x) = \frac{1}{d}\sum_{S \subseteq [d] \setminus \{j\}} \binom{d-1}{|S|}^{-1} [v_x(S \cup \{j\}) - v_x(S)]$$2

6. Theoretical Properties and Performance Guarantees

• Unbiasedness: Both constituent Monte Carlo estimators are unbiased, ensuring the final corrected estimator is unbiased for the true Shapley value of $$\phi_j(x) = \frac{1}{d}\sum_{S \subseteq [d] \setminus \{j\}} \binom{d-1}{|S|}^{-1} [v_x(S \cup \{j\}) - v_x(S)]$$3.
• Variance Reduction: For oracle $$\phi_j(x) = \frac{1}{d}\sum_{S \subseteq [d] \setminus \{j\}} \binom{d-1}{|S|}^{-1} [v_x(S \cup \{j\}) - v_x(S)]$$4, variance reduces by a factor $$\phi_j(x) = \frac{1}{d}\sum_{S \subseteq [d] \setminus \{j\}} \binom{d-1}{|S|}^{-1} [v_x(S \cup \{j\}) - v_x(S)]$$5 relative to the baseline, where $$\phi_j(x) = \frac{1}{d}\sum_{S \subseteq [d] \setminus \{j\}} \binom{d-1}{|S|}^{-1} [v_x(S \cup \{j\}) - v_x(S)]$$6 is the empirical correlation between model and surrogate Shapley estimates. In empirical settings, $$\phi_j(x) = \frac{1}{d}\sum_{S \subseteq [d] \setminus \{j\}} \binom{d-1}{|S|}^{-1} [v_x(S \cup \{j\}) - v_x(S)]$$7 is estimated from data.

7. Empirical Evaluation and Practical Considerations

ControlSHAP was evaluated on a synthetic 10-dimensional Gaussian block-correlated dataset and four UCI datasets (Bank Marketing, German Credit, Adult/Census, BRCA subtyping, $$\phi_j(x) = \frac{1}{d}\sum_{S \subseteq [d] \setminus \{j\}} \binom{d-1}{|S|}^{-1} [v_x(S \cup \{j\}) - v_x(S)]$$8–$$\phi_j(x) = \frac{1}{d}\sum_{S \subseteq [d] \setminus \{j\}} \binom{d-1}{|S|}^{-1} [v_x(S \cup \{j\}) - v_x(S)]$$9):

• Models: logistic regression, two-layer neural network, random forest.
• Protocol: For each of 40 held-out $2^d$0, 50 independent runs with $2^d$1 coalitions per run.
• Variance Reduction: Median variance reduction of 80–95% for top features in logistic regression and neural nets; up to 90% in synthetic data. For tree-based models (using finite-difference surrogates), variance reduced by 50–60%.
• Ranking Consistency: The average number of pairwise rank changes among features decreased by 30–50%.
• Bias: No detectable bias; sums of corrected Shapley values matched $2^d$2 as well or better than the basic estimators.

Practical considerations include:

• The base sample size $2^d$3 remains the main lever for controlling absolute variance. ControlSHAP effectively multiplies this by $2^d$4.
• For non-differentiable models, derivatives required for Taylor surrogates can be approximated using finite-differences.
• Correlated features require one-time $2^d$5 precomputation, but per-point computation is then constant time.
• Variance and covariance may be estimated directly, via regression for KernelSHAP, or via bootstrap.

Primary limitations are tied to the fidelity of the Taylor surrogate: for highly non-smooth $2^d$6, correlation $2^d$7 may be modest, thus diminishing variance reduction. The correlated-features surrogate assumes elliptical or Gaussian data structure.

8. Summary and Significance

ControlSHAP integrates seamlessly with existing Shapley value estimation pipelines, introducing negligible computational overhead and leveraging shared random subsets. By incorporating a locally analytic surrogate with closed-form Shapley values as a control variate, it preserves unbiasedness while affording substantial variance reductions (80–90% in practice), with broad applicability across model types and data distributions (Goldwasser et al., 2023).