--- title: "5.2.5 Feature importance: MDI, MDA, SFI and SHAP" description: "Four importance measures on a panel with known true and pure-noise features: MDI biases, MDA with purging, SFI, SHAP and groups." url: "https://pythonparatrading.com/en/python-trading-course/level-5/feature-importance-mdi-mda-sfi-shap/" lang: "en" image: "https://pythonparatrading.com/og/en/python-trading-course/level-5/feature-importance-mdi-mda-sfi-shap.png" --- # 5.2.5 · Feature importance: MDI, MDA, SFI and SHAP > Measuring which features a model uses is the most useful research tool in financial ML and one of the easiest to misread. In the level's panel it is known which features have an effect, which only inherit it from another and which are pure noise, so each method can be judged against the truth: with a noise-free target and with the real label. - Course: [Level 5: Machine learning and deep learning (Keras/PyTorch)](https://pythonparatrading.com/en/python-trading-course/level-5/) - Module: 5.2 · Feature engineering - Duration: 70 min - Prerequisites: [5.1.2 · Supervised, unsupervised and reinforcement learning; scikit-learn pipelines](https://pythonparatrading.com/en/python-trading-course/level-5/supervised-learning-scikit-learn-pipelines/), [5.2.4 · Leakage-free normalization and cross-sectional rank](https://pythonparatrading.com/en/python-trading-course/level-5/leakage-free-normalization-cross-sectional-rank/), [2.1.8 · The multiple testing problem: Bonferroni, FDR and the “factor zoo”](https://pythonparatrading.com/en/python-trading-course/level-2/multiple-testing-bonferroni-fdr-factor-zoo/) ## By the end of this lesson - Compute MDI, MDA with purged k-fold, SFI and SHAP on the same forest and explain which question each one answers. - Recognize the biases of MDI (in-sample, cardinality, sharing between correlated features) and why with a noisy label it spreads importance almost evenly. - Measure the substitution effect in MDA and SFI and correct it with importance by groups of correlated features (clustered feature importance). - Read importance with its standard error and against known noise controls, and not mistake it for proof that the feature predicts. - Check with SHAP that the explanation is of the model and not of the world: a forest that did not learn an interaction does not show it. ## What importance is for (and what it is not for) A backtest says whether a strategy made money in a sample; it does not say why. Feature importance tries to answer the second question: which variables a model uses and how much it loses without them. López de Prado (2018, chapter 8) sums it up in a working rule: the backtest is not a research tool, feature importance is. You research with importance (which phenomenon lies behind, which features are superfluous, which hypothesis deserves a new experiment) and use the backtest only once, at the end, to validate (lessons 2.5.4 and 5.4.4). There are four families of measures, and each one answers a different question: | Measure | Question | Data | Where it fails | |---|---|---|---| | MDI (*mean decrease impurity*) | How much does it reduce impurity in the tree's splits? | Training | Noise, cardinality, substitution | | MDA (*mean decrease accuracy*, permutation) | How much does the out-of-sample prediction worsen without it? | Test | Substitution, large error | | SFI (*single feature importance*) | How much does it predict on its own? | Test | Interactions, borrowed information | | SHAP | How much does it contribute to each specific prediction? | Any | Explains the model, not the data | ```mermaid flowchart TB F[Features] --> M[Model] M --> I[Importance] I --> H[Hypothesis] H --> F I -.-> B[Final backtest] ``` The problem with judging these measures on real data is that you never know what the answer is. The level 5 panel does have it. ## The test bench: known truth and two targets The panel has 400 fictitious stocks, 729 labeled weeks and 22 weekly features (without `revision_analistas`, the leak from lesson 5.2.4). The level's README publishes what each one does in the generating model: | Group | Features | Role in the generator | |---|---|---| | Own effect | `ret_1s`, `mom_12_1`, `valor_bp`, `sue`, `vol_63d`, `sentimiento` | Enter the alpha | | Only in interaction | `volumen_rel` (with `ret_1s`), `dias_desde_informe` (with `sue`) | Modulate another effect | | Inherit an effect | `ruido_4` (0.6 × momentum + noise); `ey`, `roe`, `rsi_14`, `ret_1m`, `ret_3m`, `dist_max_52s`, `vol_21d`, `log_capitalizacion` | Correlated with the true ones, no effect of their own | | Noise | `ruido_1`, `ruido_2`, `ruido_3`, `ruido_5` | None | To those 22, a control feature built in the cell itself is added here, `ruido_binario`: a coin toss, with only two values. The features enter as each week's normal rank (lesson 5.2.4) and there are two targets: - **The real label**: next week's return minus the weekly mean, in basis points. - **The true alpha** minus its mean: the predictable part of the label, without the noise. No real model has this target; it serves to see what the methods would say if the signal were not buried. ```python import numpy as np import pandas as pd from scipy import stats f = pd.read_parquet("panel5_features.parquet") e = pd.read_parquet("panel5_etiquetas.parquet") v = pd.read_parquet("panel5_verdad.parquet") d = (f.merge(e[["fecha", "ticker", "ret_sig_semana", "fin_sig_semana"]], on=["fecha", "ticker"]) .merge(v[["fecha", "ticker", "alfa_verdadero"]], on=["fecha", "ticker"]) .dropna(subset=["ret_sig_semana"]) .reset_index(drop=True)) FEATS = [c for c in f.columns if c not in ("fecha", "ticker", "sector", "revision_analistas")] sem = d.groupby("fecha") def gauss(col): n = sem[col].transform("count") return stats.norm.ppf((sem[col].rank() - 0.5) / n) X = pd.DataFrame({c: gauss(c) for c in FEATS}).fillna(0.0) rng = np.random.default_rng(0) X["ruido_binario"] = rng.integers(0, 2, len(X)).astype(float) cols = list(X.columns) def menos_media(col): # in basis points return ((d[col] - sem[col].transform("mean")) * 1e4).to_numpy() y = menos_media("ret_sig_semana") # the real label a = menos_media("alfa_verdadero") # the noise-free target print(X.shape) print(f"label std {y.std():.0f} bp · " f"alpha std {a.std():.1f} bp") ``` ```output (291600, 22) label std 366 bp · alpha std 12.7 bp ``` The label has a standard deviation of 366 bp and the alpha, of 12.7: the noise is 29 times the signal. That ratio decides almost everything that follows. The model is a random forest of 50 trees with López de Prado's choices for measuring importance: `max_features=1` (at each split one feature is chosen at random, so that none hides its correlated peers, the *masking*), large leaves (`min_samples_leaf=1000`) and 20% subsamples per tree. The out-of-sample score is the weekly IC (lesson 5.1.3), and the validation, a purged k-fold of five blocks of weeks with an embargo of four (lesson 5.4.2): with a one-week label the purge barely removes rows; the embargo prevents persistent features from connecting the end of one block with the start of the next. ```python from sklearn.ensemble import RandomForestRegressor fechas = np.sort(d["fecha"].unique()) semana = np.searchsorted(fechas, d["fecha"].to_numpy()) fin = d["fin_sig_semana"].to_numpy() inicio = d["fecha"].to_numpy() def pliegues(k=5, embargo=4): for b in np.array_split(np.arange(len(fechas)), k): te = np.isin(semana, b) ini_te, fin_te = fechas[b[0]], fin[te].max() solapa = (fin > ini_te) & (inicio <= fin_te) # purge emb = (semana > b[-1]) & (semana <= b[-1] + embargo) yield ~te & ~solapa & ~emb, te def ic_semanal(pred, obj, sem_idx): r = pd.DataFrame({"p": pred, "y": obj}).groupby(sem_idx).rank() r = r - r.groupby(sem_idx).transform("mean") num = (r["p"] * r["y"]).groupby(sem_idx).sum() den = np.sqrt((r["p"] ** 2).groupby(sem_idx).sum() * (r["y"] ** 2).groupby(sem_idx).sum()) return (num / den).mean() def bosque(n=50, muestra=0.2): return RandomForestRegressor( n_estimators=n, max_features=1, min_samples_leaf=1000, max_samples=muestra, n_jobs=4, random_state=0) for i, (tr, te) in enumerate(pliegues()): print(i, f"test {te.sum():6d} train {tr.sum():6d} " f"out {(~tr & ~te).sum():5d}") ``` ```output 0 test 58400 train 231600 out 1600 1 test 58400 train 231600 out 1600 2 test 58400 train 231600 out 1600 3 test 58400 train 231600 out 1600 4 test 58000 train 233600 out 0 ``` ## MDI: in-sample importance Each split of a tree reduces the impurity (in regression, the squared error) of its node. A feature's MDI is the sum of those reductions over all the splits that use it, weighted by the fraction of samples that pass through the node, averaged over the trees and normalized so that the total is 100%. It is free (it comes out of training) and it is what `feature_importances_` returns in scikit-learn. ```python mdi = pd.DataFrame({ "alpha": bosque().fit(X, a).feature_importances_, "label": bosque().fit(X, y).feature_importances_, }, index=cols) * 100 mdi.sort_values("alpha", ascending=False).round(2) ``` ```output alpha label ret_1s 21.60 6.47 vol_63d 12.75 6.64 valor_bp 8.17 4.76 ey 8.12 5.59 rsi_14 7.98 5.30 sue 7.60 4.63 mom_12_1 6.11 5.47 sentimiento 5.66 3.05 vol_21d 5.46 5.41 ret_1m 4.78 4.84 roe 3.98 5.06 dist_max_52s 3.06 4.36 ruido_4 1.90 4.13 log_capitalizacion 0.92 4.22 ret_3m 0.86 4.57 volumen_rel 0.53 3.87 dias_desde_informe 0.13 4.85 ruido_3 0.10 3.67 ruido_1 0.08 4.15 ruido_2 0.08 3.63 ruido_5 0.08 4.22 ruido_binario 0.02 1.13 ``` With the noise-free target, MDI resembles the truth: `ret_1s` (reversal, the strongest effect) takes 21.6%, the pure noise features 0.1% or less, and the coin 0.02%. But it already shows its first bias: it shares the effect among correlated features. `ey` (8.1%) and `rsi_14` (8.0%), with no effect of their own, beat `mom_12_1` (6.1%), which does have one, and `ret_1m` (4.8%) comes close to it: with `max_features=1`, when the tree cannot pick the true one, it picks its neighbor. With the real label, MDI separates nothing. All continuous features sit between 3.6% and 6.6%, around the 4.5% of a uniform split; `ruido_1` (4.2%) and `ruido_5` (4.2%) are within tenths of `sue` (4.6%). The reason is that MDI is measured **in sample**: if 99.9% of the label's variance is noise, the splits that reduce it most are the ones that fit noise, and any continuous feature will do. The two exceptions are `sentimiento` (3.1%), which has 32% of empty values filled with 0, and the coin (1.1%): with fewer distinct values there are fewer possible splits. It is MDI's **cardinality bias**, documented in scikit-learn's own help: it favors features with many values, even if they are noise. ## MDA: what is lost out of sample Permutation importance (MDA) is measured on data the model has not seen: it is trained on each fold, scored on the test set, one feature is shuffled in the test set (breaking its relationship with the label while keeping its distribution) and it is scored again. The drop in score is the importance; it is averaged over the folds and its standard error comes from the dispersion between them. Here the drop is measured in thousandths of IC. ```python def mda(obj, grupos=None, n=50): grupos = grupos or [[c] for c in cols] rng = np.random.default_rng(0) caidas, base = {" + ".join(g): [] for g in grupos}, [] for tr, te in pliegues(): m = bosque(n).fit(X.values[tr], obj[tr]) Xte = X.values[te] s0 = ic_semanal(m.predict(Xte), obj[te], semana[te]) base.append(s0) for g in grupos: Xp = Xte.copy() j = [cols.index(c) for c in g] Xp[:, j] = Xp[rng.permutation(len(Xp))][:, j] s = ic_semanal(m.predict(Xp), obj[te], semana[te]) caidas[" + ".join(g)].append(s0 - s) c = pd.DataFrame(caidas).T * 1000 # thousandths of IC return pd.DataFrame({"mda": c.mean(axis=1), "se": c.std(axis=1) / np.sqrt(c.shape[1])} ), np.mean(base) mda_a, base_a = mda(a) mda_y, base_y = mda(y) print(f"out-of-sample IC: alpha {base_a:.3f} · " f"label {base_y:.4f}") tabla = pd.DataFrame({"mda alpha": mda_a["mda"], "mda label": mda_y["mda"], "se label": mda_y["se"]}) tabla.sort_values("mda label", ascending=False).round(2) ``` ```output out-of-sample IC: alpha 0.652 · label 0.0148 mda alpha mda label se label ret_1s 156.82 4.04 1.09 valor_bp 25.71 3.27 0.89 ey 22.81 1.84 0.54 mom_12_1 26.92 1.67 0.58 ruido_1 0.01 1.20 0.50 ruido_binario -0.08 1.02 0.66 sentimiento 34.83 0.82 0.47 ret_1m -0.11 0.80 0.66 sue 39.69 0.70 0.60 dias_desde_informe -0.11 0.58 0.59 ruido_4 5.47 0.47 0.24 roe 2.17 0.29 0.50 dist_max_52s 8.70 0.23 0.23 volumen_rel 2.28 0.22 0.21 rsi_14 14.20 0.19 0.81 ruido_2 0.18 0.09 0.52 vol_63d 57.46 -0.04 0.49 log_capitalizacion 1.37 -0.08 0.41 ruido_5 0.29 -0.23 0.82 vol_21d 9.30 -0.41 1.02 ruido_3 0.14 -0.41 0.30 ret_3m -1.19 -0.74 0.54 ``` Out of sample, the forest achieves an IC of 0.652 with the noise-free target and of 0.0148 with the real label; the oracle's ceiling with the real label is 0.029 (lesson 5.2.4). With the noise-free target, MDA is clean: all noise features stay at zero and those with their own effect at the top. But it has the opposite bias to MDI's: `ret_1m` (−0.1) and `roe` (2.2) fall to almost nothing and `rsi_14` (14.2) sits well below `sue` or `mom_12_1`, because permuting a feature when its neighbor is still in the model barely costs anything: the model leans on the other one. It is the **substitution effect**. With the real label, MDA separates somewhat: `ret_1s` (4.0 ± 1.1), `valor_bp` (3.3 ± 0.9), `ey` (1.8 ± 0.5) and `mom_12_1` (1.7 ± 0.6) stay at the top. But the list continues with `ruido_1` (1.2 ± 0.5, $t \approx 2.4$) and `ruido_binario` (1.0 ± 0.7): two noise features ahead of `sue` and `sentimiento`, which do have an effect. And `vol_63d`, with a planted effect, comes out at zero: the label only subtracts the weekly mean, and in the panel volatile stocks have more beta, so the low-volatility premium is masked by beta (the level's README measures it). The figure at the top puts the two real-label columns side by side. With 22 features, some noise feature passing $t = 2$ is to be expected (lesson 2.1.8). That is why it pays to always put several **known noise controls** in the model, like `ruido_binario` here: their importance gives the yardstick of what chance produces, and a real feature has to sit clearly above all of them. ## SFI: each feature separately SFI removes substitution by construction: it trains a model with **a single** feature and measures its out-of-sample score, with the same validation. It is expensive (one model per feature and fold), so here each one has 10 trees and 5% subsamples. ```python def sfi(obj, n=10): res = {} for j, c in enumerate(cols): s = [ic_semanal(bosque(n, 0.05) .fit(X.values[tr][:, [j]], obj[tr]) .predict(X.values[te][:, [j]]), obj[te], semana[te]) for tr, te in pliegues()] res[c] = {"sfi": np.mean(s) * 1000, "se": np.std(s, ddof=1) / np.sqrt(len(s)) * 1000} return pd.DataFrame(res).T sfi_a, sfi_y = sfi(a), sfi(y) pd.DataFrame({"sfi alpha": sfi_a["sfi"], "sfi label": sfi_y["sfi"], "se label": sfi_y["se"]} ).sort_values("sfi alpha", ascending=False).round(1) ``` ```output sfi alpha sfi label se label ret_1s 463.8 12.5 3.0 vol_63d 300.3 2.7 1.9 sue 283.5 2.1 1.4 ey 276.4 1.9 0.9 mom_12_1 269.2 5.2 1.4 valor_bp 245.8 3.8 2.6 vol_21d 232.0 1.1 2.3 rsi_14 223.6 3.5 1.0 sentimiento 211.6 1.4 0.8 roe 185.9 2.7 1.0 ret_1m 183.1 3.3 2.2 dist_max_52s 173.9 2.6 1.4 ruido_4 159.4 -2.5 3.4 ret_3m 63.2 0.4 2.1 log_capitalizacion 51.2 0.8 1.8 ruido_3 8.7 -1.0 2.1 volumen_rel 5.5 0.7 1.7 dias_desde_informe 1.2 0.2 1.3 ruido_binario 1.1 -2.5 1.3 ruido_1 1.0 0.3 1.9 ruido_5 -1.2 1.5 1.7 ruido_2 -2.0 -1.3 2.1 ``` With the noise-free target, SFI shows the other face of substitution: `ret_1m` (183), `rsi_14` (224), `ey` (276) and `roe` (186) look almost as predictive as the true ones, because on their own they carry their neighbor's effect. And it does not see interactions: `volumen_rel` (5.5) and `dias_desde_informe` (1.2) predict nothing on their own, even though the full model uses them to modulate reversal and SUE. With the real label, SFI singles out `ret_1s` (12.5 ± 3.0) and, at a distance, `mom_12_1` (5.2 ± 1.4) and `rsi_14` (3.5 ± 1.0), which on its own carries part of the reversal. Features with no effect of their own such as `roe` (2.7 ± 1.0) sit at the level of `valor_bp` (3.8 ± 2.6) and `sue` (2.1 ± 1.4), and `ruido_4` and `ruido_binario` come out negative (−2.5): a single-feature model with no effect learns noise from the training set and pays for it in the test. ## Importance by groups of correlated features López de Prado's solution (2020, chapter 6) to substitution is to stop asking about each feature and ask about **groups**: cluster the features by their correlation and measure the MDA by permuting the whole group, so that the model cannot take refuge in a neighbor. The usual distance is $\sqrt{(1 - \rho)/2}$, which is 0 for two identical features and 1 for two opposite ones; with average linkage and a cut at 0.5 the groups come out: ```python from scipy.cluster.hierarchy import fcluster, linkage from scipy.spatial.distance import squareform C = X[FEATS].corr() dist = np.sqrt(0.5 * (1 - C.to_numpy())) Z = linkage(squareform(dist, checks=False), method="average") etiq = fcluster(Z, t=0.5, criterion="distance") grupos = {} for c, k in zip(FEATS, etiq): grupos.setdefault(k, []).append(c) grupos = list(grupos.values()) + [["ruido_binario"]] for g in grupos: if len(g) > 1: print(" + ".join(g)) print(len(grupos), "groups") ``` ```output ret_1s + ret_1m + rsi_14 ret_3m + dist_max_52s mom_12_1 + ruido_4 vol_21d + vol_63d valor_bp + ey + roe 15 groups ``` The five groups with more than one feature make economic sense: short-term reversal (`ret_1s`, `ret_1m`, `rsi_14`), medium-term return (`ret_3m`, `dist_max_52s`), momentum (`mom_12_1` and its noisy copy `ruido_4`), volatility (`vol_21d`, `vol_63d`) and value-quality (`valor_bp`, `ey`, `roe`). López de Prado proposes choosing the number of groups with a clustering algorithm using the silhouette as the criterion (ONC); a fixed threshold on a dendrogram, as here, is the simple version and is enough when the groups are this clear. ```python cmda_a, _ = mda(a, grupos) cmda_y, _ = mda(y, grupos) pd.DataFrame({"alpha": cmda_a["mda"], "label": cmda_y["mda"], "se label": cmda_y["se"]} ).sort_values("label", ascending=False).round(2) ``` ```output alpha label se label ret_1s + ret_1m + rsi_14 258.38 6.31 2.27 valor_bp + ey + roe 124.07 4.45 1.90 vol_21d + vol_63d 112.37 2.44 1.29 mom_12_1 + ruido_4 41.51 2.40 0.71 ruido_binario -0.01 1.27 0.82 log_capitalizacion 1.29 0.81 0.46 ruido_1 0.10 0.68 0.43 sue 40.23 0.58 0.63 sentimiento 34.80 0.53 0.31 dias_desde_informe 0.15 0.50 0.75 volumen_rel 2.28 0.46 0.20 ret_3m + dist_max_52s 8.43 0.28 0.66 ruido_5 0.20 -0.09 0.71 ruido_3 0.10 -0.18 0.20 ruido_2 0.08 -0.19 0.23 ``` *Fig. — MDA by groups of correlated features with the panel's real label (purged k-fold, ±2 standard errors): the four groups with a planted effect come out on top.* Points with intervals, one per group. Reversal (ret_1s, ret_1m, rsi_14), highlighted, 6.3; value-quality 4.5; volatility 2.4; momentum 2.4. Then ruido_binario 1.3, log market cap 0.8, ruido_1 0.7, sue 0.6, sentiment 0.5, days since report 0.5, volumen_rel 0.5, medium term 0.3 and ruido_5, ruido_3 and ruido_2 around −0.2. Only the reversal, value-quality and momentum intervals sit clearly above zero. With the noise-free target, the reversal group (258) collects what the individual MDA shared out badly, and the value-quality group (124) no longer depends on which of the three features the model looks at. With the real label, the four groups with a planted effect are the top four: reversal (6.3 ± 2.3), value-quality (4.5 ± 1.9), volatility (2.4 ± 1.3) and momentum (2.4 ± 0.7). Behind them, the coin (1.3 ± 0.8) is a reminder of the size of the measurement's noise. ## SHAP: the contribution to each prediction Shapley values come from cooperative game theory: they split the “payoff” of a coalition among its players according to their average contribution over all possible coalitions. SHAP (Lundberg and Lee, 2017) applies them to a prediction: the difference between a row's prediction and the mean is split among the features, and the sum of the values is exactly that difference. For trees, `shap.TreeExplainer` computes them exactly and fast. The usual global importance is the mean of the absolute values: ```python import shap tr, te = list(pliegues())[-1] m_y = bosque().fit(X.values[tr], y[tr]) muestra = np.random.default_rng(1).choice(np.flatnonzero(te), 3000, replace=False) sv = shap.TreeExplainer(m_y).shap_values(X.values[muestra]) shap_y = pd.Series(np.abs(sv).mean(axis=0), index=cols) shap_y.sort_values(ascending=False).round(2).head(10) ``` ```output ret_1s 1.02 rsi_14 0.94 ret_1m 0.83 ruido_1 0.71 mom_12_1 0.67 ruido_5 0.66 valor_bp 0.64 ey 0.59 roe 0.53 vol_63d 0.51 dtype: float64 ``` Computed on test rows of the last fold with the forest trained on the real label, the mean SHAP ranks `ret_1s` first, but puts `ruido_1` (0.71) and `ruido_5` (0.66) in the top 10 ahead of `valor_bp`: SHAP explains what the model does, and this model, trained on noise, uses noise. It is not a measure of how much a feature predicts, but of how much it moves the predictions. SHAP's advantage is in the detail: it says how a feature's contribution depends on its value and on those of the others. The generator plants an interaction: reversal is multiplied by $(1 - 0.5\cdot g(\text{volumen\_rel}))$, stronger when volume is low. For a term $c\cdot x\cdot(1 - 0.5z)$, the Shapley value of $x$ is $c\cdot x\cdot(1 - 0.25z)$: half of the cross term is attributed to $z$. The slope of the SHAP of `ret_1s` in each tercile of `volumen_rel`, relative to the middle tercile, can be compared with that exact answer: ```python import lightgbm as lgb m_a = bosque().fit(X.values[tr], a[tr]) m_l = lgb.LGBMRegressor( n_estimators=300, learning_rate=0.05, num_leaves=15, min_child_samples=500, subsample=0.5, subsample_freq=1, random_state=0, n_jobs=4, deterministic=True, force_row_wise=True, verbose=-1).fit(X.values[tr], a[tr]) xs = X.iloc[muestra] j = cols.index("ret_1s") tercio = pd.cut(xs["volumen_rel"], [-np.inf, -0.43, 0.43, np.inf], labels=["low", "medium", "high"]) def pendientes(modelo): sv = shap.TreeExplainer(modelo).shap_values(X.values[muestra]) p = pd.Series({k: np.polyfit(xs["ret_1s"][tercio == k], sv[(tercio == k).to_numpy(), j], 1)[0] for k in ["low", "medium", "high"]}) return p / p["medium"] media_g = xs.groupby(tercio, observed=True)["volumen_rel"].mean() t = pd.DataFrame({"forest": pendientes(m_a), "LightGBM": pendientes(m_l), "exact SHAP": 1 - 0.25 * media_g, "true effect": 1 - 0.5 * media_g}) t.round(2) ``` ```output forest LightGBM exact SHAP true effect low 1.05 1.25 1.27 1.55 medium 1.00 1.00 1.00 1.00 high 0.91 0.72 0.72 0.45 ``` *Fig. — Slope of the ret_1s SHAP value in each volumen_rel tercile, relative to the middle tercile, with the noise-free target: LightGBM, the max_features=1 forest, the exact Shapley value and the generator's true effect.* Four lines from low to high relative volume. True effect: 1.55, 1 and 0.45. Exact Shapley, which assigns half of the interaction to volumen_rel: 1.27, 1 and 0.72. LightGBM, highlighted: 1.25, 1 and 0.72, almost on top of the exact value. Forest: 1.05, 1 and 0.91, nearly flat: it barely learned the interaction. The LightGBM (trees of up to 15 leaves, which can combine features along the same path) gives 1.25 and 0.72, right on the exact Shapley value (1.27 and 0.72). The forest of the previous sections gives 1.05 and 0.91: with one random feature per split and leaves of 1,000 rows, it barely learned the interaction, and SHAP reflects that faithfully. The explanation is of the model, not of the data: if the model did not learn something, SHAP is not going to show it, and if it learned noise, it will show it as if it were signal. ## How to use all this ```mermaid flowchart TB A[Candidates] --> G[Groups] G --> P[MDA with purging] P --> C{Above controls?} C -- no --> D[Discard] C -- yes --> S[SHAP: shape] S --> V[Final validation] ``` - **Never MDI alone.** With financial labels it is almost uniform and biased by cardinality. If it is used, it should be with the target of a problem that has signal or as a diagnostic of the model, not of the data. - **MDA with the model's same honest validation** (purging, embargo) and with the economic score (IC, Sharpe of the signal), not with raw accuracy. - **Groups before individual features** when there are strong correlations, and always with noise controls. - **Standard error and multiple testing.** An MDA with $t = 2$ among 22 features is not a finding. - **Selection overfits too.** Choosing features by their importance over the whole sample and then validating on the same sample is a selection bias: the choice is made within each training window or a period that nobody has looked at is set aside (lessons 5.4.3 and 5.4.4). ## Exercises 1. Repeat the MDI with the noise-free target and `max_features=None` (all features at each split). What happens to the features correlated with `ret_1s`? **Solution.** ```python todas = RandomForestRegressor( n_estimators=20, max_features=None, min_samples_leaf=1000, max_samples=0.2, n_jobs=4, random_state=0).fit(X, a) comp = pd.DataFrame({ "max_features=1": mdi["alpha"], "max_features=None": todas.feature_importances_ * 100}) comp.loc[["ret_1s", "ret_1m", "rsi_14", "valor_bp", "ey", "roe", "mom_12_1", "ruido_4"]].round(2) ``` ```output max_features=1 max_features=None ret_1s 21.60 52.35 ret_1m 4.78 0.00 rsi_14 7.98 0.00 valor_bp 8.17 2.03 ey 8.12 8.60 roe 3.98 0.00 mom_12_1 6.11 1.84 ruido_4 1.90 0.00 ``` With all features available at each split, the tree always picks the one that reduces impurity most and the neighbors disappear: `ret_1s` goes from 21.6% to 52.4% and `ret_1m`, `rsi_14`, `roe` and `ruido_4` fall to zero. It is *masking*. And the true one does not always win: in the value-quality group `ey` keeps almost everything (8.6%) and `valor_bp`, the one with the planted effect, drops to 2.0%, because at each split the one that best fits that sample wins, not the causal one. With `max_features=1` the split between correlated features is fairer, at the cost of all of them looking somewhat important; that is why group importance is the underlying answer. 2. Add `revision_analistas` (the leak from lesson 5.2.4) to the model and compute the MDA of four features with the real label. Does importance detect the leak? **Solution.** ```python X["revision_analistas"] = gauss("revision_analistas") X["revision_analistas"] = X["revision_analistas"].fillna(0.0) cols = list(X.columns) mda_fuga, base_fuga = mda(y, [["revision_analistas"], ["ret_1s"], ["valor_bp"], ["ruido_1"]]) print(f"out-of-sample IC with the leak: {base_fuga:.3f}") mda_fuga.round(2) ``` ```output out-of-sample IC with the leak: 0.157 mda se revision_analistas 155.61 2.89 ret_1s 0.68 0.53 valor_bp 0.53 0.24 ruido_1 -0.15 0.13 ``` The leak comes out a hundred times more important than any other feature (156 thousandths of IC against less than 1) and the forest's out-of-sample IC rises from 0.015 to 0.157, five times the oracle's ceiling. Importance does not detect the leak; it points to it, for those who know how to read it. A feature that takes almost all the importance in a problem with so little signal is not a finding: it is the first place to look for a date error. With the leak in the model, moreover, the rest of the features lose importance (that of `ret_1s` drops from 4.0 to 0.7): the model explains with the leak part of what it used to explain with them. 3. You have 300 candidate features for a weekly stock selection model and ten years of data. Describe a protocol for keeping a few of them without fooling yourself. **Solution.** First, a date audit of each feature (lesson 5.2.4): no importance measure detects a leak. Then, transform into ranks within each date and cluster by correlation (with distance $\sqrt{(1-\rho)/2}$ and a criterion such as the silhouette) to work with tens of groups and not with hundreds of features. Add several noise controls (random features with the same persistence as the real ones: one of i.i.d. noise and a slow AR(1) one). Measure the MDA by groups with purged k-fold and the economic score, only on the first seven years; keep the groups clearly above the best control, with the number of tests in mind (lesson 2.1.8). Within each group, choose the simplest feature or the one with the best SFI. Look with SHAP at the shape of each effect in case it contradicts the economic hypothesis. And validate the final model only once on the last three years, which were not used to choose anything, with a Deflated Sharpe over all the attempts (lesson 5.4.4). ## Self-assessment 1. With the panel's real label, a forest's MDI gives `ruido_1` (i.i.d. noise) 4.2% of the importance and `sue` (planted effect) 4.6%. Why does MDI not separate them? - A. Because `sue` has no effect - B. Because MDI is measured in sample: with a signal-to-noise ratio of 1 to 30, the trees use any continuous feature to fit noise and the split comes out almost uniform - C. Because the features are normal ranks - D. Because there are not enough trees *Answer B: Because MDI is measured in sample: with a signal-to-noise ratio of 1 to 30, the trees use any continuous feature to fit noise and the split comes out almost uniform.* MDI adds up the impurity reduction of the splits on the training data. If almost all the variation in the label is noise, the splits that reduce it most are the ones that fit noise, and any continuous feature serves equally well. With the noise-free target (the true alpha) the MDI of `ruido_1` is 0.08%. 2. `ruido_binario` (a coin toss) receives 1.1% of the MDI with the real label; `ruido_1`, 4.2%. Both are noise. Where does the difference come from? - A. From the coin having a negative effect - B. From the generator's seed - C. From the forest using `max_features=1` - D. From MDI's cardinality bias: a feature with two values offers a single possible split and a continuous one, thousands, so the continuous one has more chances to fit noise *Answer D: From MDI's cardinality bias: a feature with two values offers a single possible split and a continuous one, thousands, so the continuous one has more chances to fit noise.* It is MDI's documented bias in favor of features with many distinct values. With the noise-free target, both have an MDI of almost zero. 3. With the noise-free target, `ret_1m` (no effect of its own, correlation 0.48 with `ret_1s`) has an MDA of −0.1 and an SFI of 183. What does the difference indicate? - A. That the SFI is miscalculated - B. That `ret_1m` has a nonlinear effect - C. A substitution effect: alone, `ret_1m` carries information borrowed from `ret_1s`; with `ret_1s` in the model, permuting it costs nothing because the model uses the original - D. That MDA needs more permutations *Answer C: A substitution effect: alone, `ret_1m` carries information borrowed from `ret_1s`; with `ret_1s` in the model, permuting it costs nothing because the model uses the original.* MDA measures what the model loses without the feature given the rest; SFI, what the feature predicts on its own. Both readings are correct and answer different questions; group importance brings `ret_1s`, `ret_1m` and `rsi_14` together and assigns them the joint effect. 4. With the real label, `ruido_1` has an MDA of 1.20 ± 0.50 thousandths of IC (t ≈ 2.4) in purged k-fold. What do you conclude? - A. That with 22 features and such a noisy label, some noise feature will come out “significant” by chance: you have to compare with known noise controls and correct for multiple testing - B. That `ruido_1` has a small effect and should be kept - C. That the purged k-fold has a leak - D. That MDA is always biased upward *Answer A: That with 22 features and such a noisy label, some noise feature will come out “significant” by chance: you have to compare with known noise controls and correct for multiple testing.* It is a textbook false positive. With several known noise features in the model, their distribution of importances gives the yardstick: a real feature has to sit clearly above all of them, and not just pass t = 2 (lesson 2.1.8). 5. A LightGBM trained on the noise-free target gives `ret_1s` relative SHAP slopes of 1.25 with low volume and 0.72 with high volume; the true effect is 1.55 and 0.45. Why do they not match? - A. Because SHAP splits the interaction between the two features involved: half goes to `volumen_rel`, and the expected values for `ret_1s` are 1.27 and 0.72, which is what it measures - B. Because the model did not learn the interaction - C. Because SHAP does not work with trees - D. Because the sample is small *Answer A: Because SHAP splits the interaction between the two features involved: half goes to `volumen_rel`, and the expected values for `ret_1s` are 1.27 and 0.72, which is what it measures.* For a term c·x·(1 − 0.5·z), the Shapley value of x is c·x·(1 − 0.25·z): the other half of the cross term is attributed to z. The forest with `max_features=1` and large leaves, on the other hand, barely learned the interaction (1.05 and 0.91), and SHAP reflects it: it explains the model, not the data. ## Key points - MDI, MDA, SFI and SHAP answer different questions: how much the model uses each feature in training (MDI), how much it loses out of sample without it (MDA), how much it predicts on its own (SFI) and how much it contributes to each prediction (SHAP). - With the panel's real label, MDI spreads importance almost evenly between true features and noise (from 3.6% to 6.6% against the uniform 4.5%) and penalizes the binary feature: it is in-sample and biased by cardinality. - MDA with purged k-fold and the IC as the score separates the true features somewhat (`ret_1s` 4.0 ± 1.1 thousandths of IC), but with a large standard error: pure noise comes out with t ≈ 2.4. Importance is read with its error and against noise controls. - Correlated features substitute for each other: MDA undervalues them and SFI lends them their neighbor's effect. Group importance (correlation, average linkage, threshold 0.5) assigns the effect to the group. - SHAP explains the model: a LightGBM that learned the reversal × volume interaction shows it with the exact Shapley slopes; a forest that did not learn it does not. - No importance measure detects a leak: `revision_analistas` comes out a hundred times more important than the rest. A disproportionate importance is a reason for an audit, not for celebration. ## Resources ### Documentation - [Lundberg and Lee (2017): A Unified Approach to Interpreting Model Predictions (SHAP), arXiv:1705.07874](https://arxiv.org/abs/1705.07874) · arxiv.org - [scikit-learn: permutation importance vs MDI in a random forest](https://scikit-learn.org/stable/auto_examples/inspection/plot_permutation_importance.html) · scikit-learn.org - [SHAP: TreeExplainer](https://shap.readthedocs.io/en/latest/generated/shap.TreeExplainer.html) · shap.readthedocs.io ### Books - [López de Prado (2018): Advances in Financial Machine Learning, Wiley (ch. 8, feature importance)](https://www.wiley.com/en-us/Advances+in+Financial+Machine+Learning-p-9781119482086) · wiley.com - [López de Prado (2020): Machine Learning for Asset Managers, Cambridge Elements (ch. 6, clustered feature importance)](https://www.cambridge.org/core/elements/machine-learning-for-asset-managers/6D9211305EA2E425D33A9F38D0AE3545) · cambridge.org ### Data - [Level 5 data: panel features (Parquet)](https://pythonparatrading.com/curso/datos/nivel-5/panel5_features.parquet) · Download · panel5_features.parquet - [Level 5 data: panel true alpha (Parquet)](https://pythonparatrading.com/curso/datos/nivel-5/panel5_verdad.parquet) · Download · panel5_verdad.parquet Previous lesson: [5.2.4 · Leakage-free normalization and cross-sectional rank](https://pythonparatrading.com/en/python-trading-course/level-5/leakage-free-normalization-cross-sectional-rank/) Next lesson: [5.3.1 · Fixed-horizon labeling and its problems](https://pythonparatrading.com/en/python-trading-course/level-5/fixed-horizon-labeling/)