Hooker Mentch Zhou 2021 - Partial Dependence Alternatives and Conclusions

Summary

Fixing effect plots is harder than fixing importance scores. §6 discusses three options.

  1. Restrict ICE curves to the range of supported by the data. This needs an estimated conditional distribution.
  2. Generalized functional ANOVA (Hooker 2007): find closest to in under the joint feature density. This needs and a high-dimensional .
  3. Additive-model distillation: fit a GAM to or to the data (Lou et al. 2013; Tan et al. 2018). This gives a coherent, extrapolation-free picture.

Simply averaging only the “supported” parts of ICE curves can be strongly misleading. §7 concludes that every sound alternative costs additional computation or modelling. Condition-and-relearn is best for importance. Local methods (LIME, saliency) avoid extrapolation but give no global summary, and counterfactual explanations may share the problem.

Overview

§5 fixed importance by staying on the data distribution. For plots of effects this is “less straightforward” (p. 16), because a PD curve is defined by varying over its whole range at every observation.

Main Content

Option 1: Restricted ICE plots (§6, p. 16)

  • In Fig. 3, each ICE line is drawn solid only over the high-density range of given the other features. This gives local information and also shows how strongly depends on the other covariates.
  • That range must be determined somehow. Goldstein et al. (2015) suggested marking the observations on ICE plots, but this does not show over what range the curve rests on nearby data.
  • The paper used ±2 standard deviations within the Gaussian copula that generated the data. Without that knowledge, the bounds require an estimated conditional distribution.

Caution: Averaging only the supported ICE segments (§6, p. 16)

Fig. 6 might suggest modifying PD to average only the solid parts of each ICE line. The authors warn this “could be strongly misleading”, because changes due to get attributed to . Different observations contribute at different values, and their values also differ.

Option 2: Generalized functional ANOVA (§6, pp. 16-18)

Definition: Distribution-weighted additive projection (Hooker 2007) (§6, p. 16)

Choose functions and to minimize

where approximates the joint feature distribution and is an unknown function of all features except . The fitted plays the role of a partial dependence function, but it is weighted by where the data actually are. It is the generalized version of the functional ANOVA in eq. (1) (^def-fanova).

  • Hooker (2007) minimized a quasi-Monte Carlo approximation of the integral, which needs an estimate of . That estimate is unlikely to be accurate in high dimensions.
  • Chastaing et al. (2012) added smoothing methods (generalized Hoeffding-Sobol decomposition).
  • Mentch & Hooker (2017) used similar structures to test for feature interactions in random forests.
  • Key obstacle: and the high-dimensional must be estimated jointly, which needs specialized learning methods.

Option 3: Additive models / distillation (§6, p. 18)

  • Lou et al. (2013) and Tan et al. (2018) fit a generalized additive model (Wood 2006) containing only the one-dimensional terms on the right of eq. (1).
  • The GAM can be fit to the data directly, or to the values of as model distillation.
  • The linear-model fits in Fig. 3 belong to this class and provide “a coherent, extrapolation-free, representation of the underlying function.”

Definition: Additive model of the one-dimensional terms of eq. (1) (§6, p. 18; formal statement is the vault's rendering)

The source describes this option in words: a GAM “that includes only one-dimensional functions in the right hand side of (1)“. Written out, the functional ANOVA (1) is , and the additive model keeps only the first two parts:

with each a smooth univariate function (e.g. a penalized spline, Wood 2006), centred for identifiability, e.g. . There are two ways to fit it:

  • To the data: .
  • As model distillation of a black box : . The black box’s predictions at the observed are the response.

The effect display for feature is the curve .

Algorithm: Additive-model distillation of a fitted model (§6, p. 18; steps are the vault's rendering)

  1. Fit the flexible model (forest, network, boosting) as usual.
  2. Compute at the observed feature vectors only, . Optionally use held-out or extra in-distribution points. Never use points on a grid or from permuted columns.
  3. Fit to the pairs by penalized least squares.
  4. Plot each over the observed range of . Report the fidelity of the distillation, e.g. of against on the observed points, so readers know how much of (interactions) the display leaves out.

Why this avoids extrapolation, and how it relates to Option 2.

  • The fitting criterion only involves (or ) at observed . Unlike PD, nothing is evaluated at off-support combinations such as (large , small ) in the copula design. This is the sense of “coherent, extrapolation-free” (p. 18).
  • Comparing it with the distribution-weighted additive projection above, the additive model is the generalized functional ANOVA criterion with two simplifications: the integral against is replaced by the empirical average over the sample, which is automatically weighted by where the data lie, and the high-dimensional is restricted to the additive . That restriction is what removes the “key problem” of jointly estimating a high-dimensional (p. 18). The price is that interactions are not shown.
  • Fig. 3’s linear fits are the special case . Their “effect curves” are the Theorem 1 lines of slope .
  • Vault gloss, not in the paper: the same compensation effect applies. With highly correlated the backfitting splits their shared signal between and , and these components are less stable than those of independent features. Additive displays avoid extrapolation but not collinearity. Compare ^def-compensation.

Conclusions (§7, p. 18)

  1. Naïve permutation methods can mislead even in simple examples, and extrapolation by flexible models is a significant source of error. The precise bias depends on the learner, the feature dependence and the response function; Bénard et al. (2021) analyse it for random forests (Sobol-MDA).
  2. All alternatives need more computation.
    • still permutes but refits. Mentch & Hooker (2016) found that keeping the feature dimension gives a better statistical comparison than dropping and refitting ().
    • Condition-and-relearn is better still, because the replacements keep at least approximate dependence on the remaining features.
    • Berrett et al. (2020), the conditional permutation test, is a hybrid: it permutes, but samples permutations non-uniformly according to the conditional distribution.
    • Non-refitting approaches need a model for ; Liu & Zheng (2018) build auto-encoding knockoff generators.
  3. Generalized Sobol indices (Hooker 2007) and Shapley values (Owen 2014; Lundberg & Lee 2017) are further options. If they are not computed carefully, they can show the same extrapolation bias.
  4. PD and ICE plots raise the same concern. Local methods (LIME, saliency) avoid extrapolation but give no global picture, for example of a threshold function. Counterfactual explanations (Wachter et al. 2017) may extrapolate in a similar way; the authors did not study this.
  5. Additive-model displays (Tan et al. 2018) represented overall function behaviour better than combinations of local explanations. The diagnostics of Hooker (2007) still need specialized methods.

Limitations and scope (as stated or implied in the paper)

  • The paper studies bias only. Variability and confidence intervals are left to other work, such as U-statistics (Ishwaran & Lu 2019) and model class reliance (§1).
  • The analytic conditional distribution in the simulations is a convenience that practice rarely has (§5).
  • Theorems 2-3 are exact for squared error and linear models respectively. Beyond these, the measures may differ (§5).
  • The extrapolation explanations are given for random forests and single-hidden-layer neural networks only.

Examples

Example: response curves from a flexible MMM surrogate (generic illustration)

A team fits gradient-boosted trees to 104 weeks of sales and plots the partial dependence of sales on TV spend to read off saturation. TV spend never exceeded its 90th percentile in weeks with low search spend. The top of the PD curve therefore averages predictions for (high TV, low search) combinations the model never saw, and its flattening or rising may reflect how trees extrapolate rather than saturation. Two safer options follow from §6. One is to show ICE only over the supported range of TV given the other inputs. The other is to fit an explicitly additive saturation model and read curves from that; see Shape (Saturation) Effects.

Connections

See Also