LUCID: Local Understanding of Covariate Importance with Conditional Density Estimation

Yixiao Wang; Srikar Katta; Jon Donnelly; Emanuele Borgonovo; Cynthia Rudin

cs.LGOct 26, 2026

LUCID: Local Understanding of Covariate Importance with Conditional Density Estimation

Authors: Yixiao Wang, Srikar Katta, Jon Donnelly, Emanuele Borgonovo, Cynthia Rudin

Organizations: Duke University, Durham, USA · Johns Hopkins University, Baltimore, USA · Bocconi University, Milan, Italy

Abstract

Computing local-level variable importance (VI) is critical for understanding machine learning models. Local permutation-based VI metrics offer a simple yet powerful approach, measuring how much a model's output changes when perturbing its inputs. However, existing perturbation strategies often ignore feature dependence or account for it only through restricted approximations. Perturbing a feature without fully respecting this dependence can generate unrealistic counterfactuals, thereby evaluating the model in regions unsupported by the observed data distribution. We address this gap by introducing a new model-agnostic local VI metric and estimator, Local Understanding of Covariate Importance with conditional Density estimation (LUCID). LUCID perturbs features of interest using plausible values drawn from the full conditional distribution given the remaining observed features. We establish theoretical connections between LUCID and existing VI measures, and provide an estimator with convergence guarantees and uncertainty quantification. Finally, we demonstrate LUCID's stable performance and practical utility across simulated and real-world examples.

Figures & tables

Figure 1: Two motivations for LUCID on the Bike Sharing dataset ( Fanaee-T and Gama, 2014 ) . (a) LUCID values for temperature and humidity LUCID (solid lines) vary substantially across hours, while their corresponding global VI values (dashed lines) obscure this heterogeneity. (b) For temperature , each split violin compares the estimated conditional replacement distribution from DR-BART ( Orlandi et al., 2021 , left;) with the marginal replacement distribution (right). Conditional densities adapt across seasons, whereas the marginal distribution does not.

Figure 2: Three examples of conditional distributions of perturbed outputs

g(Z^{\prime},\bm{x})

. The purple dashed line marks the mean

\mathbb{E}\!\left[g(Z^{\prime},\bm{x})\right]

, the yellow arrow shows the square root of the displacement , while the purple arrow illustrates the conditional standard deviation, the square root of the dispersion term.

Category	Method
Attribution	SHAP ( Lundberg and Lee, 2017 )
	LIME ( Ribeiro et al., 2016 )
	Cond. SHAP ( Aas et al., 2021 ; Jullum et al., 2025 )
Removal	local version of LOCO ( Lei et al., 2018 )
Permutation	ICI ( Casalicchio et al., 2019 )
	CMR ( Fisher et al., 2019 )

Table 1: Methods considered in our experiments.

	DR-BART	ARF	cs-PFI	CMR	GCMR	CLIQUE	Permutation
DGP 1	0.022	0.144	0.008	0.590	0.543	0.001	0.008
DGP 2	0.004	0.005	0.015	0.011	0.011	4.460	4.431
DGP 3	0.449	0.454	1.173	2.685	2.485	2.060	2.065
DGP 4	0.401	1.298	83.329	1.809	4.315	83.673	83.135
DGP 5	0.066	0.185	0.036	0.849	0.851	0.251	0.269
Mean rank	1.800	2.800	3.400	5.200	5.400	4.600	4.800

Table 2: MSE between estimated and analytical LUCID, averaged over all evaluated instances of

X_{1}

. The lowest MSE (from unrounded MSE values) in each DGP is highlighted in purple bold . Mean rank averages within-DGP rankings; lower is better.

Figure 3: Discontinuity. Local importance for $X_{1}$ across methods. Red boxes mark representative nonzero importance in the inactive region. Solid and dashed curves denote analytical LUCID and CPFI, respectively.

Figure 4: Heteroskedastic. Resampling distributions of

X_{1}

at five

X_{2}

quantile points across samplers.

Figure 5: LUCID estimation MSE as independent irrelevant features are added across five DGPs.

	DGP 1	DGP 2	DGP 3	DGP 4	DGP 5
90% $C$	98.43%	91.70%	91.34%	93.20%	97.12%
90% $W$	0.766	0.155	0.944	0.969	0.956
95% $C$	99.37%	95.76%	92.94%	95.64%	98.52%
95% $W$	0.910	0.184	1.123	1.152	1.138

Table 3: Empirical coverage (

C

) and width (

W

) of LUCID posterior credible intervals for active

X_{1}

observations across five DGPs.

Figure 6: Local importance and conditional replacements for hour. (a) LUCID for hour varies substantially over the day and between working and non-working days (solid lines). (b–d) Split violins compare the estimated conditional distribution of replacement hour values using DR-BART (left) with the marginal permutation distribution (right).

Figure 7: weather-group LUCID (a) The purple line shows the joint LUCID of temperature, humidity, and season, while the orange dashed line shows the sum of their individual LUCIDs; the shaded region highlights their gap. (b) Group gaps are shown for the

(T,H,S)

(T,S)

(T,H)

, and

(H,S)

groups.

Appendix figures & tables35 assets

Supplementary material from the paper’s appendix.

Appendix

Notation	Definition
$Z,\ z$	feature of interest and its observed value.
$\bm{X}=(X_{1},\ldots,X_{k}),\ \bm{x}$	Remaining covariates and their observed values; $k$ is the number of conditioning covariates.
$\bm{W}=(Z,\bm{X}),\ \bm{w}=(z,\bm{x})$	Complete feature vector and one observation.
$Z^{\prime},\ z^{\prime}$	conditional replacement and its integration variable.
$\bm{Z},\ \bm{z},\ \bm{Z}^{\prime},\ \bm{z}^{\prime}$	Feature group, its observed vector, its joint conditional replacement, and the vector integration variable.
$P_{Z,\bm{X}},\ P_{\bm{X}},\ P_{Z}$	True joint covariate distribution and its marginals.

Appendix

Table 4: Notation for LUCID and its theoretical guarantees.

DGP	Data-generating process	Analytical LUCID	CPFI
Discontinuity	$X_{1},X_{2},X_{3}\stackrel{{\scriptstyle\mathrm{iid}}}{{\sim}}\mathcal{N}(0,1)$ , $g(\bm{x})=x_{1}\mathbf{1}_{x_{3}>0}+x_{2}\mathbf{1}_{x_{3}\leq 0}$ .	$I_{1}=(x_{1}^{2}+1)\mathbf{1}_{x_{3}>0}$ ; $I_{2}=(x_{2}^{2}+1)\mathbf{1}_{x_{3}\leq 0}$ ; $I_{3}=\frac{1}{2}(x_{1}-x_{2})^{2}$ .	$\mathrm{CPFI}_{1}=1$ ; $\mathrm{CPFI}_{2}=1$ ; $\mathrm{CPFI}_{3}=1$ .
Correlated	$(X_{1},X_{2})^{\top}\sim\mathcal{N}_{2}(0,\Sigma_{\rho})$ , $\Sigma_{\rho,11}=\Sigma_{\rho,22}=1$ , $\Sigma_{\rho,12}=\rho$ , $\rho=0.9$ ; $g(\bm{x})=x_{1}$ .	$I_{1}=1-\rho^{2}+(x_{1}-\rho x_{2})^{2}$ ; $I_{2}=0$ .	$\mathrm{CPFI}_{1}=2(1-\rho^{2})$ ; $\mathrm{CPFI}_{2}=0$ .
Heteroskedastic	$X_{2}\sim\mathcal{N}(0,1)$ , $X_{1}\mid X_{2}=x_{2}\sim\mathcal{N}(0,x_{2}^{2})$ ; $g(\bm{x})=x_{1}$ .	$I_{1}=x_{1}^{2}+x_{2}^{2}$ ; $I_{2}=0$ .	$\mathrm{CPFI}_{1}=2$ ; $\mathrm{CPFI}_{2}=0$ .
Interacted discontinuity	$X_{2},X_{3}\stackrel{{\scriptstyle\mathrm{iid}}}{{\sim}}\mathcal{N}(0,1)$ , $\mu=4(\mathbf{1}_{X_{2}>0}\oplus\mathbf{1}_{X_{3}>0})-2$ ; $X_{1}\mid X_{2},X_{3}\sim\mathcal{N}(\mu,1)$ ; $g(\bm{x})=x_{1}$ .	$I_{1}=1+(x_{1}-\mu)^{2}$ ; $I_{2}=I_{3}=0$ .	$\mathrm{CPFI}_{1}=2$ ; $\mathrm{CPFI}_{2}=\mathrm{CPFI}_{3}=0$ .
Extreme tail	$X_{2}\sim\mathcal{N}(0,1)$ , $X_{1}\mid X_{2}=x_{2}\sim\mathcal{N}(0,v(x_{2}))$ ; $v(x_{2})=1+\mathbf{1}_{x_{2}>0}$ , $H(x_{2})=\mathbf{1}_{\|x_{2}\|>\Phi^{-1}(.99)}$ ; $g(\bm{x})=x_{1}+20H(x_{2})$ .	$I_{1}=x_{1}^{2}+v(x_{2})$ ; $I_{2}=8+384H(x_{2})$ .	$\mathrm{CPFI}_{1}=3$ ; $\mathrm{CPFI}_{2}=15.68$ .

Appendix

Table 5: Synthetic DGPs and analytical LUCID and analytical CPFI under squared loss.

Method	DGP 1	DGP 2	DGP 3	DGP 4	DGP 5
DR-BART	$0.0225\pm 0.0036$	$\mathbf{0.0029\pm 0.0006}$	$\mathbf{0.3668\pm 0.1095}$	$\mathbf{0.5092\pm 0.1293}$	$0.0865\pm 0.0122$
ARF	$0.1187\pm 0.0147$	$0.0041\pm 0.0006$	$0.6582\pm 0.2888$	$1.2859\pm 0.1738$	$0.2063\pm 0.0194$
cs-PFI	$0.0064\pm 0.0019$	$0.0116\pm 0.0016$	$1.0804\pm 0.2386$	$78.0902\pm 5.7078$	$\mathbf{0.0510\pm 0.0125}$
CLIQUE	$\mathbf{0.0008\pm 0.0007}$	$4.5720\pm 0.2497$	$1.9922\pm 0.1122$	$80.3271\pm 1.7660$	$0.2536\pm 0.0027$
CMR	$0.5556\pm 0.0586$	$0.0106\pm 0.0010$	$2.5597\pm 0.2128$	$1.5994\pm 0.2780$	$0.8806\pm 0.0583$
GCMR	$0.5447\pm 0.0603$	$0.0108\pm 0.0010$	$2.4399\pm 0.1563$	$3.9210\pm 0.3313$	$0.8883\pm 0.0590$

Appendix

Table 6: MSE of estimated LUCID for

X_{1}

across 50 independent random seeds. Each entry reports mean

\pm

standard deviation across seeds. Lower values indicate more accurate and more stable estimation.

Method	Mean rank
DR-BART	1.80
ARF	2.84
cs-PFI	3.29
CLIQUE	4.60
Permutation	4.80
CMR	5.16

Appendix

Table 7: Mean MSE rank for

X_{1}

across the five DGPs and 50 independent random seeds, giving

250

paired DGP–seed comparisons. Lower is better.

Figure 8: Extreme tail. Match and mismatch regions for the analytical and fitted LUCID structures in DGP 5. Match and mismatch regions between LUCID and local LOCO feature rankings in DGP 5, using analytical and estimated importance measures.

Figure 9: Estimated conditional replacement distributions for

X_{1}

in the discontinuity DGP at five empirical quantiles of

X_{2}

, with

X_{3}=1

. The analytical conditional distribution is shown as Optimal.

Figure 10: Estimated conditional replacement distributions for

X_{1}

in the correlated-features DGP at five empirical quantiles of

X_{2}

Figure 11: Estimated conditional replacement distributions for

X_{1}

in the heteroskedastic DGP at five empirical quantiles of

X_{2}

Figure 12: Estimated conditional replacement distributions for

X_{1}

in the interacted-discontinuity DGP at five empirical quantiles of

X_{2}

, with

X_{3}=1

Figure 13: Estimated conditional replacement distributions for

X_{1}

in the extreme-tail DGP at five empirical quantiles of

X_{2}

Figure 14: Complete pointwise LUCID results for the regime-dependent DGP. Estimated local variable importance values from all conditional sampling methods are compared with the corresponding analytical LUCID across evaluation observations.

Figure 15: Complete pointwise LUCID results for the strong conditional dependence DGP. Estimated local variable importance values from all conditional sampling methods are compared with the corresponding analytical LUCID across evaluation observations.

Figure 16: Complete pointwise LUCID results for the heteroskedastic conditional dependence DGP. Estimated local variable importance values from all conditional sampling methods are compared with the corresponding analytical LUCID across evaluation observations.

Figure 17: Complete pointwise LUCID results for the interaction conditional dependence DGP. Estimated local variable importance values from all conditional sampling methods are compared with the corresponding analytical LUCID across evaluation observations.

Figure 18: Complete pointwise LUCID results for the extreme-tail DGP. Estimated local variable importance values from all conditional sampling methods are compared with the corresponding analytical LUCID across evaluation observations.

	50%		80%		90%		95%
DGP	$C$	$W$	$C$	$W$	$C$	$W$	$C$	$W$
DGP 1	0.686	0.311	0.948	0.596	0.984	0.766	0.994	0.910
DGP 2	0.519	0.063	0.828	0.121	0.917	0.155	0.958	0.184
DGP 3	0.640	0.382	0.871	0.734	0.913	0.944	0.929	1.123
DGP 4	0.606	0.393	0.872	0.753	0.932	0.969	0.956	1.152
DGP 5	0.676	0.387	0.929	0.742	0.971	0.956	0.985	1.138

Appendix

Table 8: Empirical coverage (

C

) and average width (

W

) of posterior LUCID credible intervals for observations where

X_{1}

is active in the prediction function across nominal levels.

Figure 19: 50% posterior credible intervals. Pointwise posterior LUCID credible intervals for $X_{1}$ across the five DGPs at the $50\%$ nominal level. Posterior means are compared with analytical LUCID values, with covered and missed analytical values indicated separately.

Figure 20: 80% posterior credible intervals. Pointwise posterior LUCID credible intervals for $X_{1}$ across the five DGPs at the $80\%$ nominal level. Posterior means are compared with analytical LUCID values, with covered and missed analytical values indicated separately.

Figure 21: 90% posterior credible intervals. Pointwise posterior LUCID credible intervals for $X_{1}$ across the five DGPs at the $90\%$ nominal level. Posterior means are compared with analytical LUCID values, with covered and missed analytical values indicated separately.

Figure 22: 95% posterior credible intervals. Pointwise posterior LUCID credible intervals for $X_{1}$ across the five DGPs at the $95\%$ nominal level. Posterior means are compared with analytical LUCID values, with covered and missed analytical values indicated separately.

Setting	Conditional distribution of $X_{1}$	Analytical LUCID	CPFI
Independent	$h(\bm{X}_{-1})=0$ , $\epsilon_{x}\sim\mathcal{N}(0,1)$	$(x_{2}+1)^{2}(1+x_{1}^{2})$	$4$


Linear	$h(\bm{X}_{-1})=X_{2}$ , $\epsilon_{x}\sim\mathcal{N}(0,1)$	$(x_{2}+1)^{2}\{1+(x_{1}-x_{2})^{2}\}$	$4$

Appendix

Table 9: Data-generating processes from the cs-PFI simulation and the corresponding analytical LUCID and CPFI for

X_{1}

under squared loss.

Independent
Setting	DR-BART	ARF	cs-PFI	CMR	GCMR	CLIQUE	Perm.
$n=300,\ p=10$	3.30	13.23	1.29	21.94	15.68	0.13	0.95
$n=300,\ p=90$	5.03	4.04	0.31	17.31	10.32	0.26	0.25
$n=3000,\ p=10$	1.21	15.65	0.22	10.89	16.05	0.014	0.21
$n=3000,\ p=90$	1.21	7.44	0.10	9.53	8.23	0.004	0.10
Linear

Appendix

Table 10: MSE between estimated and analytical LUCID for

X_{1}

in the cs-PFI simulation settings. Lower values indicate more accurate estimation. The lowest MSE in each setting is shown in bold.

Figure 23: Complete individual LUCID results for the Bike Sharing data. For temperature, hour, working day, humidity, and season, panel (a) shows mean LUCID across hour of day separately for working and non-working days. Panels (b)–(d) compare the DR-BART conditional replacement distributions across season, working-day status, and hour periods with the corresponding marginal permutation distributions.

Figure 24: Feature relationships and individual LUCID in the Bike Sharing data. The left panel summarizes pairwise correlations and each feature’s predictability $R_{j}^{2}$ from the remaining covariates. The right panel shows observation-level LUCID values, colored by the observed feature value. Temperature is highly predictable from the remaining features, while hour and working day have very low predictability, providing additional context for interpreting their individual LUCID values.

Group	$\bm{R_{G}^{2}}$	$\bm{R_{G,\mathrm{within}}^{2}}$
Season + temperature	0.198	0.450
Temperature + humidity	0.482	0.190
Humidity + season	0.298	0.054
Season + humidity + temperature	0.095	0.410

Appendix

Table 11: Group outside predictability

R_{G}^{2}

and within-group predictability

R_{G,\mathrm{within}}^{2}

for candidate weather-related feature groups in the Bike Sharing data.

Figure 25: Dispersion–displacement decomposition of hour LUCID. Observation-level hour LUCID is decomposed into the dispersion term $V$ and the displacement term $D$ , the percentages indicate their respective shares of the mean LUCID across observations. Mean LUCID composition: $V=51.3\%$ , $D=48.7\%$ of the mean LUCID. Color indicates which component contributes more strongly, while the solid and dashed curves summarize working and non-working days, respectively.

Figure 26: Conditional replacement distributions for hour. Replacement distributions produced by DR-BART, ARF, cs-PFI, CMR, GCMR, and marginal permutation across season, working-day status, and hour periods. Hour is weakly predictable from the remaining features, and its conditional replacement distribution is correspondingly broad and relatively stable across most contexts.

Figure 27: Conditional replacement distributions for humidity. Replacement distributions across season, working-day status, and hour periods. Several conditional methods capture systematic changes in humidity across covariate contexts, particularly across hour periods, whereas marginal permutation uses the same unconditional distribution throughout. CMR also exhibits broader support than the naturally bounded humidity feature.

Figure 28: Conditional replacement distributions for temperature. Replacement distributions across season, working-day status, and hour periods. The strongest distributional shifts occur across seasons, where conditional samplers adapt to the changing temperature distribution while marginal permutation remains fixed across contexts.

Figure 29: Conditional replacement distributions for working day. Replacement distributions across season, working-day status, and hour periods. Most conditional methods preserve the binary structure of the feature, while CMR can place mass outside the natural $\{0,1\}$ support. The conditional distribution varies relatively little across the displayed contexts.

Figure 30: Conditional replacement distributions for season. Replacement distributions across observed season, working-day status, and hour periods. Conditional methods preserve differing degrees of context dependence in the categorical replacement distribution, whereas marginal permutation does not condition on the observed context. CMR again produces values outside the natural categorical support.

Dataset	Samples ( $n$ )	Features ( $p$ )
Heart	3,658	15
German	1,000	60
HELOC	9,871	23
GMSC	102,209	10

Appendix

Table 12: OpenXAI benchmark datasets used in the scalability experiments.

Method	Heart	German	HELOC	GMSC
Permutation	$8.81\pm 1.16$	$6.47\pm 0.64$	$15.67\pm 2.21$	$87.18\pm 21.91$
CMR	$5.89\pm 0.03$	$5.35\pm 0.28$	$14.99\pm 3.37$	$92.92\pm 23.30$
GCMR	$7.37\pm 1.85$	$6.84\pm 0.14$	$16.33\pm 2.77$	$106.18\pm 24.79$
local LOCO	$5.41\pm 0.88$	$4.97\pm 0.15$	$5.60\pm 1.29$	$4.70\pm 0.18$
CLIQUE	$7.65\pm 1.68$	$5.50\pm 0.34$	$14.61\pm 2.87$	$117.08\pm 3.76$
ARF	$43.84\pm 5.88$	$30.10\pm 3.69$	$287.20\pm 25.11$	$2451.36\pm 816.10$

Appendix

Table 13: Average runtime per feature-specific task in seconds (mean

\pm

standard deviation across five folds). Each task includes conditional-sampler fitting and LUCID computation.

Method	Heart	German	HELOC	GMSC
Permutation	$2.20\pm 0.29$	$6.47\pm 0.64$	$6.01\pm 0.85$	$14.53\pm 3.65$
CMR	$1.47\pm 0.01$	$5.35\pm 0.28$	$5.75\pm 1.29$	$15.49\pm 3.88$
GCMR	$1.84\pm 0.46$	$6.84\pm 0.14$	$6.26\pm 1.06$	$17.70\pm 4.13$
local LOCO	$1.35\pm 0.22$	$4.97\pm 0.15$	$2.15\pm 0.50$	$0.78\pm 0.03$
CLIQUE	$1.91\pm 0.42$	$5.50\pm 0.34$	$5.60\pm 1.10$	$19.51\pm 0.63$
ARF	$10.96\pm 1.47$	$30.10\pm 3.69$	$110.09\pm 9.63$	$408.56\pm 136.02$

Appendix

Table 14: Serial-equivalent runtime per fold in minutes (mean

\pm

standard deviation across five folds), obtained by summing all feature-specific jobs. This represents the runtime required when feature computations are not parallelized.

Method	Heart	German	HELOC	GMSC
Permutation	$0.338\pm 0.001$	$0.329\pm 0.000$	$0.352\pm 0.000$	$0.365\pm 0.004$
CMR	$0.335\pm 0.000$	$0.327\pm 0.002$	$0.354\pm 0.003$	$0.371\pm 0.002$
GCMR	$0.339\pm 0.003$	$0.332\pm 0.001$	$0.358\pm 0.001$	$0.375\pm 0.005$
local LOCO	$0.320\pm 0.003$	$0.320\pm 0.002$	$0.325\pm 0.003$	$0.341\pm 0.000$
CLIQUE	$0.336\pm 0.003$	$0.325\pm 0.002$	$0.351\pm 0.004$	$0.367\pm 0.001$
ARF	$0.996\pm 0.006$	$2.072\pm 0.030$	$1.705\pm 0.004$	$0.985\pm 0.010$

Appendix

Table 15: Peak memory usage in GiB (mean

\pm

standard deviation across five folds). For each fold, the maximum peak memory among feature-specific jobs is reported.

Perturbation-Based Feature Attribution

Perturbation-based feature attribution explains model predictions by altering inputs and measuring how those changes affect outputs.

See all 41 papers

CID: Measuring Feature Importance Through Counterfactual Distributions
Nov 19, 2025Eddie Conti, Álvaro Parafita, Axel BrandoFeature Attribution Explainability Evaluation
TAIES Group, HPES Lab, Barcelona Supercomputing Center (BSC)
Regional Explanations: Bridging Local and Global Variable Importance
Apr 13, 2026Salim I. Amoukou, Nicolas J-B. BrunelFeature Attribution Explainable Artificial Intelligence
J.P. Morgan AI Research · LaMME, ENSIIE, University Paris Saclay & Capgemini Invent
One Permutation Is All You Need: Fast, Deterministic Feature Importance and Model Stress-Testing
Dec 15, 2025Albert DoradorPerturbation-Based Feature Attribution Algorithmic Auditing

Nov 19, 2025J/K move · Enter open · S save

CID: Measuring Feature Importance Through Counterfactual Distributions

Eddie Conti, Álvaro Parafita, Axel Brando

TAIES Group, HPES Lab, Barcelona Supercomputing Center (BSC)

Feature Attribution Explainability Evaluation

Open paper AskPDF arXiv

Assessing the importance of individual features in Machine Learning is critical to understand the model's decision-making process. While numerous methods exist, the lack of a definitive ground truth for comparison highlights the need for alternative, well-founded measures. This paper introduces a novel post-hoc local feature importance method called Counterfactual Importance Distribution (CID). We generate two sets of positive and negative counterfactuals, model their distributions using Kernel Density Estimation, and rank features based on a distributional dissimilarity measure. This measure, grounded in a rigorous mathematical framework, satisfies key properties required to function as a valid metric. We showcase the effectiveness of our method by comparing with well-established local feature importance explainers. Our method not only offers complementary perspectives to existing approaches, but also improves performance on faithfulness metrics (both for comprehensiveness and sufficiency), resulting in more faithful explanations of the system. These results highlight its potential as a valuable tool for model analysis. Link to repository: https://github.com/EddieConti/CID

Abstract

Perturbation-Based Feature Attribution

Explore similar work

CID: Measuring Feature Importance Through Counterfactual Distributions

Regional Explanations: Bridging Local and Global Variable Importance

One Permutation Is All You Need: Fast, Deterministic Feature Importance and Model Stress-Testing