Explainability

SHAP, LIME, permutation importance, per-atom contribution maps, fragment analysis, partial dependence and counterfactuals.

For a QSAR model the interesting question is not “which feature index mattered” but “which part of this molecule mattered”, which is what the atom-level tools answer.

Permutation importance

Model-agnostic, and the only importance measure here that reflects the model’s actual predictive reliance rather than its internal structure:

>>> from qsarkit.explainability import PermutationImportance
>>> from qsarkit.models import QSARRegressor
>>> X, y = demo_fingerprints(256), DEMO_Y
>>> model = QSARRegressor("rf", random_state=0).fit(X, y)
>>> importance = PermutationImportance(n_repeats=5, random_state=0).fit(model, X, y)
>>> importance.importances_mean_.shape
(256,)

Measure it on held-out data. Permuting a feature on the training set reports how much the model memorized through it, which is not the same quantity and is usually larger.

Per-atom contributions

The chemist-facing explanation: mask each atom in turn, re-predict, and attribute the change to that atom.

>>> from qsarkit.explainability import AtomicContributionMap
>>> from qsarkit.representation import MorganFingerprint
>>> explainer = AtomicContributionMap(model, MorganFingerprint(n_bits=256))
>>> result = explainer.explain(demo_mols[0])
>>> result.weights.shape
(9,)

One weight per heavy atom, ready for RDKit’s similarity maps.

Note

Masking one atom of a symmetric ring changes little, because the remaining atoms set the same bits. That is not a defect of the method — it is a true statement about the model: no single one of those atoms is necessary, because the others are redundant with it. Read flat weights across a ring as “this ring matters as a unit”.

SHAP and LIME

Both need the explainability extra:

from qsarkit.explainability import LIMEExplainer, SHAPExplainer

shap_values = SHAPExplainer(model).explain(X)      # needs qsarkit-learn[explainability]

SHAP’s TreeExplainer is exact and fast for tree ensembles; KernelExplainer is model-agnostic and slow enough that you will want to subsample. LIME fits a local surrogate, which makes it cheap and local-only — a LIME explanation says nothing about a different molecule.

Fragments and partial dependence

>>> from qsarkit.explainability import FragmentContributionAnalyzer, PartialDependence
>>> analyzer = FragmentContributionAnalyzer(model, MorganFingerprint(n_bits=256))
>>> hasattr(analyzer, "analyze")
True
>>> pd = PartialDependence(model, grid_resolution=10)
>>> hasattr(pd, "compute")
True

Fragment contributions aggregate atom-level attributions over chemically meaningful groups, which is usually more actionable than either the atom or the bit.

API

Model interpretation – OECD validation principle 5.

Two families, answering different questions:

For a fingerprint model the second family is the one that matters: 2048 anonymous bits are not an explanation, but “the nitro group costs you a log unit” is.

Examples

>>> from rdkit import Chem
>>> from sklearn.ensemble import RandomForestRegressor
>>> from qsarkit.explainability import AtomicContributionMap
>>> from qsarkit.representation import MorganFingerprint
>>> fp = MorganFingerprint(n_bits=64)
>>> mols = [Chem.MolFromSmiles(s) for s in ("CCO", "CCN", "c1ccccc1", "CCC")]
>>> model = RandomForestRegressor(n_estimators=5, random_state=0).fit(
...     fp.transform(mols), [1.0, 2.0, 3.0, 4.0]
... )
>>> AtomicContributionMap(model, fp).explain(mols[0]).weights.shape
(3,)

References

class qsarkit.explainability.AttributionAtomMapper(fingerprint, distribution='uniform')[source]

Bases: object

Turn a SHAP or LIME explainer’s output into atom-level weights.

Wraps any per-feature attribution and the fingerprint that produced those features, so a model explanation can be shown on the structure rather than as bit indices.

Parameters:

  • fingerprint (Any) – The transformer whose bits the attributions refer to. Its radius, n_bits and (where present) use_features attributes are read so the bit numbering matches – typically a MorganFingerprint.

  • distribution (Literal['uniform', 'center', 'radius_weighted']) – Passed to bit_weights_to_atom_weights().

Variables:

Examples

>>> import numpy as np
>>> from qsarkit.explainability import AttributionAtomMapper
>>> from qsarkit.models import QSARRegressor
>>> from qsarkit.representation import MorganFingerprint
>>> from rdkit import Chem
>>> mols = [Chem.MolFromSmiles(s) for s in
...         ("CC(=O)Nc1ccccc1", "CC(=O)Nc1ccc(Cl)cc1", "CCO", "CCN")]
>>> fingerprint = MorganFingerprint(radius=2, n_bits=256)
>>> X = fingerprint.transform(mols)
>>> model = QSARRegressor("rf", random_state=0).fit(X, [6.2, 6.7, 5.0, 5.1])
>>> mapper = AttributionAtomMapper(fingerprint)

Any per-bit vector maps onto atoms. Here we use permutation importance, which works for every backend and needs no optional dependency:

>>> from qsarkit.explainability import PermutationImportance
>>> importance = PermutationImportance(n_repeats=2, random_state=0).fit(
...     model, X, [6.2, 6.7, 5.0, 5.1])
>>> weights = mapper.atom_weights(mols[1], importance.importances_mean_)
>>> weights.shape == (mols[1].GetNumAtoms(),)
True

Symmetry-equivalent atoms receive identical weight, which is a useful correctness check on the mapping:

>>> chlorobenzene = Chem.MolFromSmiles("Clc1ccccc1")
>>> import numpy as np
>>> bits = np.ones(256)
>>> symmetric = mapper.atom_weights(chlorobenzene, bits)
>>> bool(np.isclose(symmetric[2], symmetric[6]))
True

The collision rate says how much to trust the picture:

>>> rate = mapper.collision_rate(mols[1])
>>> 0.0 <= rate <= 1.0
True

References

atom_weights(mol, bit_weights)[source]

Map one molecule’s per-bit attributions onto its atoms.

Parameters:
Return type:

ndarray[tuple[Any, ...], dtype[double]]

from_shap(mol, explainer, X, index=0)[source]

Atom weights from a SHAPExplainer.

Parameters:

  • mol (Any) – The molecule corresponding to row index of X.

  • explainer (Any) – A fitted explainer.

  • X (Union[Buffer, _SupportsArray[dtype[Any]], _NestedSequence[_SupportsArray[dtype[Any]]], complex, bytes, str, _NestedSequence[complex | bytes | str]]) – The feature matrix the explanation is computed on.

  • index (int) – Which row of X to explain.

Return type:

ndarray[tuple[Any, ...], dtype[double]]

Raises:

OptionalDependencyError – If shap is not installed.

from_lime(mol, explanation)[source]

Atom weights from a LIME explanation of one molecule.

Parameters:

  • mol (Any)

  • explanation (Any) – Either a LIME Explanation (as_map() is read), or any mapping of feature index -> weight.

Return type:

ndarray[tuple[Any, ...], dtype[double]]

collision_rate(mol)[source]

Fraction of this molecule’s set bits shared by several environments.

A bit set by more than one environment carries attribution that cannot be assigned to a single substructure, so a high rate means the atom-level picture is blurred. Widening the fingerprint lowers it.

Parameters:

mol (Any)

Returns:

In [0, 1]. 0.0 when the molecule sets no bits.

Return type:

float

qsarkit.explainability.bit_atom_environments(mol, radius=2, n_bits=2048, use_features=False)[source]

Which atom environments set each bit of a molecule’s Morgan fingerprint.

Parameters:

  • mol (Any) – The molecule.

  • radius (int) – Morgan radius. Must match the fingerprint the attributions came from, or the bit numbering will not correspond.

  • n_bits (int) – Fingerprint length. Must likewise match.

  • use_features (bool) – Use the feature-based (FCFP) invariants rather than connectivity.

Returns:

Bit index -> tuple of (central_atom_index, radius) pairs, one per environment in this molecule that sets that bit.

Return type:

Dict[int, Tuple[Tuple[int, int], ...]]

Examples

>>> from rdkit import Chem
>>> from qsarkit.explainability import bit_atom_environments
>>> mol = Chem.MolFromSmiles("CCO")
>>> environments = bit_atom_environments(mol, n_bits=256)
>>> len(environments) > 0
True
>>> all(isinstance(bit, int) for bit in environments)
True

Every reported central atom is a real atom of the molecule:

>>> centres = {atom for envs in environments.values() for atom, _ in envs}
>>> max(centres) < mol.GetNumAtoms()
True

References

qsarkit.explainability.bit_weights_to_atom_weights(mol, bit_weights, radius=2, n_bits=2048, use_features=False, distribution='uniform')[source]

Spread per-bit attributions over the atoms that produced each bit.

Parameters:

  • mol (Any) – The molecule the attributions were computed for.

  • bit_weights (Union[Buffer, _SupportsArray[dtype[Any]], _NestedSequence[_SupportsArray[dtype[Any]]], complex, bytes, str, _NestedSequence[complex | bytes | str]]) – One attribution per fingerprint bit – a row of SHAP values, a LIME coefficient vector, or any per-feature importance.

  • radius (int) – Morgan radius used to compute bit_weights.

  • n_bits (int) – Fingerprint length used to compute bit_weights.

  • use_features (bool) – Whether those were FCFP rather than ECFP bits.

  • distribution (Literal['uniform', 'center', 'radius_weighted']) –

    How a bit’s weight is divided among its environment’s atoms:

    "uniform"

    Split equally. Neutral, and the right default.

    "center"

    All of it to the central atom. Sharper pictures, but it overstates the centre of a large environment.

    "radius_weighted"

    Split equally, then divide by radius + 1, so a bit describing a tight environment carries more weight per atom than one describing a diffuse environment.

Returns:

Per-atom weight, ready for draw_atom_weights().

Return type:

ndarray[tuple[Any, ...], dtype[double]]

Raises:

ValueError – If bit_weights is not one-dimensional of length n_bits.

Examples

>>> import numpy as np
>>> from rdkit import Chem
>>> from qsarkit.explainability import (
...     bit_atom_environments, bit_weights_to_atom_weights)
>>> mol = Chem.MolFromSmiles("CC(=O)Nc1ccc(Cl)cc1")
>>> weights = np.zeros(256)
>>> environments = bit_atom_environments(mol, n_bits=256)
>>> chlorine = [a.GetIdx() for a in mol.GetAtoms() if a.GetSymbol() == "Cl"][0]
>>> # Attribute to one bit centred on the chlorine, radius 0.
>>> for bit, envs in environments.items():
...     if (chlorine, 0) in envs:
...         weights[bit] = 1.0
>>> atom_weights = bit_weights_to_atom_weights(mol, weights, n_bits=256)
>>> int(np.argmax(atom_weights)) == chlorine
True

Only the atoms in that environment receive weight:

>>> int((atom_weights != 0).sum())
1

References

qsarkit.explainability.draw_atom_weights(mol, atom_weights, size=(400, 400), fmt='svg', normalize=True, contour_lines=10)[source]

Render atom weights on the structure as an RDKit similarity map.

The standard cheminformatics depiction: a green-to-pink field over the 2D structure, positive contributions in one colour and negative in the other.

Parameters:

  • mol (Any) – The molecule. A 2D conformer is computed if it has none.

  • atom_weights (Union[Buffer, _SupportsArray[dtype[Any]], _NestedSequence[_SupportsArray[dtype[Any]]], complex, bytes, str, _NestedSequence[complex | bytes | str]]) – Per-atom weight, e.g. from AttributionAtomMapper.

  • size (Tuple[int, int]) – Image size in pixels.

  • fmt (Literal['svg', 'png']) – "svg" returns a string, which renders inline in a notebook and embeds in HTML. "png" returns bytes.

  • normalize (bool) – Scale the weights so the largest absolute value maps to the end of the colour scale. Keeps the picture readable regardless of the attribution’s units; turn it off to compare two molecules on one absolute scale.

  • contour_lines (int) – Number of contour lines drawn.

Returns:

SVG text, or PNG bytes.

Return type:

Any

Raises:

ValueError – If atom_weights does not have one entry per atom, or fmt is not recognized.

Examples

>>> import numpy as np
>>> from rdkit import Chem
>>> from qsarkit.explainability import draw_atom_weights
>>> mol = Chem.MolFromSmiles("CC(=O)Nc1ccc(Cl)cc1")
>>> weights = np.linspace(-1, 1, mol.GetNumAtoms())
>>> svg = draw_atom_weights(mol, weights)
>>> svg.lstrip().startswith("<?xml") or svg.lstrip().startswith("<svg")
True
>>> isinstance(draw_atom_weights(mol, weights, fmt="png"), bytes)
True

A weight per atom is required, so a mismatched vector is caught rather than silently misaligned:

>>> draw_atom_weights(mol, [0.1, 0.2])
Traceback (most recent call last):
    ...
ValueError: atom_weights has 2 entries but the molecule has 11 atoms.

References

class qsarkit.explainability.PermutationImportance(n_repeats=10, scoring=None, random_state=None, n_jobs=None)[source]

Bases: object

Feature importance by measuring the damage from shuffling a column.

Model-agnostic and honest: a feature matters exactly to the degree that destroying its relationship with the target degrades held-out performance. Unlike a tree model’s built-in feature_importances_, which is computed on training data and is biased toward high-cardinality features, this is measured on data the model has not seen.

The caveat worth knowing: with correlated descriptors — which is the normal situation in QSAR — permuting one of a correlated pair leaves the model able to recover the signal from its partner, so both look unimportant. Cluster correlated descriptors before interpreting, or read the result as “importance given the others are present”.

Parameters:

  • n_repeats (int) – Shuffles per feature. More repeats reduce the variance of the estimate.

  • scoring (Optional[str]) – scikit-learn scorer name. Defaults to the estimator’s own score.

  • random_state (Optional[int]) – Seed.

  • n_jobs (Optional[int]) – Parallel jobs.

Variables:

  • importances_mean (ndarray) – Mean drop in score per feature.

  • importances_std (ndarray) – Standard deviation across repeats.

Examples

>>> from sklearn.datasets import make_regression
>>> from sklearn.ensemble import RandomForestRegressor
>>> X, y = make_regression(n_samples=60, n_features=5, n_informative=2,
...                        random_state=0)
>>> model = RandomForestRegressor(n_estimators=10, random_state=0).fit(X, y)
>>> imp = PermutationImportance(random_state=0).fit(model, X, y)
>>> imp.importances_mean_.shape
(5,)

References

importances_mean_: ndarray[tuple[Any, ...], dtype[float64]]
importances_std_: ndarray[tuple[Any, ...], dtype[float64]]
fit(estimator, X, y, feature_names=None)[source]

Measure importance on the supplied (ideally held-out) data.

Parameters:
Return type:

PermutationImportance

to_dataframe(top_n=None)[source]

Importances as a ranked table.

Parameters:

top_n (Optional[int]) – Keep only the top_n most important features.

Returns:

Columns feature, importance, std, descending.

Return type:

DataFrame

plot(top_n=20)[source]

Horizontal bar chart of the most important features.

Parameters:

top_n (int)

Return type:

Figure

class qsarkit.explainability.SHAPExplainer(model, explainer_type='auto', background=None, n_background=100, random_state=None)[source]

Bases: object

SHAP values: the game-theoretic attribution of a prediction.

SHAP assigns each feature the payoff it contributes to a prediction, averaged over all orderings in which features could be added. That construction gives it the properties ad-hoc attributions lack — the contributions sum exactly to the prediction minus the base value (local accuracy), and a feature the model ignores always gets zero.

The explainer is chosen from the model type: TreeExplainer for forests and boosted trees (exact and fast), LinearExplainer for linear models, and KernelExplainer otherwise (model-agnostic but slow, so it samples the background set).

Requires the explainability extra.

Parameters:

  • model (Any)

  • explainer_type (Literal['auto', 'tree', 'linear', 'kernel']) – Which SHAP explainer to use.

  • background (Union[Buffer, _SupportsArray[dtype[Any]], _NestedSequence[_SupportsArray[dtype[Any]]], complex, bytes, str, _NestedSequence[complex | bytes | str], None]) – Background dataset for the kernel/linear explainers. A sample of the training data; 100 rows is usually enough.

  • n_background (int) – How many background rows to sample when background is a full training set.

  • random_state (Optional[int])

Examples

>>> import pytest
>>> shap = pytest.importorskip("shap")

References

property explainer: Any

The lazily-constructed SHAP explainer.

shap_values(X)[source]

SHAP values for each sample and feature.

Parameters:

X (Union[Buffer, _SupportsArray[dtype[Any]], _NestedSequence[_SupportsArray[dtype[Any]]], complex, bytes, str, _NestedSequence[complex | bytes | str]])

Returns:

For multiclass models, the values for the positive class.

Return type:

ndarray[tuple[Any, ...], dtype[double]]

global_importance(X, feature_names=None)[source]

Mean absolute SHAP value per feature — a global ranking.

Parameters:
Returns:

Columns feature, importance, descending.

Return type:

DataFrame

explain_one(x, feature_names=None, top_n=10)[source]

Per-feature contributions to a single prediction.

Parameters:
Returns:

Columns feature, value, shap_value, ordered by absolute contribution.

Return type:

DataFrame

plot_importance(X, feature_names=None, top_n=20)[source]

Bar chart of mean absolute SHAP value per feature.

Parameters:
Return type:

Figure

class qsarkit.explainability.LIMEExplainer(model, training_data, feature_names=None, mode='regression', n_samples=5000, random_state=None)[source]

Bases: object

LIME: explain one prediction with a local surrogate model.

Perturbs the molecule’s descriptors, records what the model predicts for each perturbation, and fits a sparse linear model to that local neighbourhood. The surrogate’s coefficients are the explanation.

Compared with SHAP, LIME is faster and easier to read but its explanations are not guaranteed to be self-consistent: the answer depends on the perturbation distribution and kernel width, and re-running can give a different story. Prefer SHAP where the attribution has to be defensible.

Requires the explainability extra.

Parameters:

Examples

>>> import pytest
>>> lime = pytest.importorskip("lime")

References

  • Ribeiro, M. T., Singh, S. & Guestrin, C. (2016). “Why Should I Trust You?: Explaining the Predictions of Any Classifier.” KDD 2016, 1135-1144. https://doi.org/10.1145/2939672.2939778

  • Alvarez-Melis, D. & Jaakkola, T. S. (2018). “On the Robustness of Interpretability Methods.” arXiv:1806.08049. https://arxiv.org/abs/1806.08049

property explainer: Any

The lazily-constructed LIME tabular explainer.

explain_one(x, top_n=10)[source]

Explain a single prediction.

Parameters:
Returns:

Columns feature, weight, ordered by absolute weight.

Return type:

DataFrame

class qsarkit.explainability.PartialDependence(model, grid_resolution=50)[source]

Bases: object

Partial dependence: the model’s average response to one descriptor.

Sweeps a descriptor across its range, averaging the model’s prediction over the observed distribution of the others. The resulting curve shows the shape of the model’s dependence — whether logP acts linearly, saturates, or has an optimum — which a single importance number cannot express.

Its known blind spot is extrapolation: averaging over the marginal distribution evaluates the model at descriptor combinations that never occur (a molecule with MW 100 and 40 rotatable bonds), so read the curve only across the range where the data are dense.

Parameters:

  • model (Any)

  • grid_resolution (int) – Points sampled across each feature’s range.

Examples

>>> from sklearn.datasets import make_regression
>>> from sklearn.ensemble import RandomForestRegressor
>>> X, y = make_regression(n_samples=50, n_features=4, random_state=0)
>>> model = RandomForestRegressor(n_estimators=5, random_state=0).fit(X, y)
>>> grid, avg = PartialDependence(model).compute(X, feature=0)
>>> grid.shape == avg.shape
True

References

compute(X, feature)[source]

Partial-dependence curve for one feature.

Parameters:
Return type:

tuple[ndarray[tuple[Any, ...], dtype[double]], ndarray[tuple[Any, ...], dtype[double]]]

Returns:

  • grid (ndarray of shape (grid_resolution,)) – Feature values swept.

  • average (ndarray of shape (grid_resolution,)) – Mean prediction at each grid point.

plot(X, feature, feature_name=None)[source]

Plot the partial-dependence curve.

Parameters:
Return type:

Figure

class qsarkit.explainability.AtomicContribution(weights, prediction, smiles='')[source]

Bases: object

Per-atom attribution for one molecule.

Variables:

  • weights (ndarray of shape (n_atoms,)) – Signed contribution of each atom. Positive means removing the atom would lower the prediction.

  • prediction (float) – The model’s prediction for the intact molecule.

  • smiles (str) – Canonical SMILES of the molecule explained.

weights: ndarray[tuple[Any, ...], dtype[float64]]
prediction: float
smiles: str
property most_positive: int

Index of the atom contributing most positively.

property most_negative: int

Index of the atom contributing most negatively.

class qsarkit.explainability.AtomicContributionMap(model, fingerprint, use_proba=False)[source]

Bases: object

Attribute a prediction to individual atoms by fingerprint masking.

Implements the similarity-map algorithm of Riniker and Landrum: for each atom, the fingerprint is recomputed with that atom removed from the environment, and the change in prediction is the atom’s contribution. Every bit set by an atom is switched off, so the difference measures exactly what that atom’s substructures were worth to the model.

This is what turns a fingerprint model — otherwise a black box over 2048 anonymous bits — into something a chemist can act on: the output maps directly onto the structure, showing which part of the molecule the model is actually responding to. It is the practical route to OECD principle 5 for fingerprint QSAR.

One property of the method is worth knowing: on a symmetric molecule every atom receives the same weight, and on a redundant one the weights are near zero. Masking a single carbon of benzene leaves the other five still generating aromatic bits, so the prediction barely moves and no atom looks responsible. That is the honest answer — the signal is carried by the ring as a whole, not by any one atom — but it means a flat weight vector indicates redundancy rather than irrelevance. Use FragmentContributionAnalyzer when the question is about a group rather than an atom.

Parameters:

  • model (Any) – Must expose predict, and predict_proba if use_proba=True.

  • fingerprint (Any) – Maps Iterable[Mol] -> ndarray. Any qsarkit.representation fingerprint works.

  • use_proba (bool) – Explain predict_proba(...)[:, 1] rather than predict, which gives a smoothly varying signal for classifiers instead of a step function.

Examples

>>> from rdkit import Chem
>>> from sklearn.ensemble import RandomForestRegressor
>>> from qsarkit.representation import MorganFingerprint
>>> fp = MorganFingerprint(n_bits=64)
>>> mols = [Chem.MolFromSmiles(s) for s in ("CCO", "CCN", "c1ccccc1", "CCC")]
>>> model = RandomForestRegressor(n_estimators=5, random_state=0).fit(
...     fp.transform(mols), [1.0, 2.0, 3.0, 4.0]
... )
>>> result = AtomicContributionMap(model, fp).explain(mols[0])
>>> result.weights.shape
(3,)

References

explain(mol)[source]

Compute per-atom contributions for one molecule.

Parameters:

mol (Any)

Return type:

AtomicContribution

Raises:

ValueError – If mol is None.

transform(mols)[source]

Explain a batch of molecules.

Parameters:

mols (Sequence[Any])

Return type:

List[AtomicContribution]

to_similarity_map_weights(mol)[source]

Weights in the form RDKit’s SimilarityMaps drawing expects.

Parameters:

mol (Any)

Returns:

One weight per atom, ready for SimilarityMaps.GetSimilarityMapFromWeights.

Return type:

List[float]

plot(mol)[source]

Bar chart of per-atom contributions, labelled by element.

Parameters:

mol (Any)

Return type:

Figure

class qsarkit.explainability.FragmentContributionAnalyzer(model, fingerprint, fragments=None, use_proba=False)[source]

Bases: object

Attribute predictions to chemically meaningful fragments.

Per-atom weights answer “which atoms matter”; chemists think in groups. This aggregates atomic contributions over substructures — a supplied SMARTS list, or the molecule’s own BRICS fragments — so the output is “the nitro group contributes -1.2 log units” rather than a list of atom indices.

Parameters:

  • model (Any)

  • fingerprint (Any)

  • fragments (Optional[Dict[str, str]]) – name -> SMARTS to attribute against. When None, BRICS decomposition is used to find the molecule’s own fragments.

  • use_proba (bool)

Examples

>>> from rdkit import Chem
>>> from sklearn.ensemble import RandomForestRegressor
>>> from qsarkit.representation import MorganFingerprint
>>> fp = MorganFingerprint(n_bits=64)
>>> mols = [Chem.MolFromSmiles(s) for s in ("CCO", "CCN", "c1ccccc1O", "CCC")]
>>> model = RandomForestRegressor(n_estimators=5, random_state=0).fit(
...     fp.transform(mols), [1.0, 2.0, 3.0, 4.0]
... )
>>> analyzer = FragmentContributionAnalyzer(
...     model, fp, fragments={"hydroxyl": "[OX2H]"}
... )
>>> df = analyzer.analyze(mols[0])
>>> "hydroxyl" in list(df["fragment"])
True

References

analyze(mol)[source]

Contribution of each fragment to the prediction.

Parameters:

mol (Any)

Returns:

Columns fragment, n_atoms, contribution (summed over the fragment’s atoms), mean_contribution, sorted by descending absolute contribution. Fragments not present in the molecule are omitted.

Return type:

DataFrame

class qsarkit.explainability.CounterfactualExplainer(model, fingerprint, delta=1.0, min_similarity=0.4, use_proba=False)[source]

Bases: object

Explain a prediction by finding the nearest molecule the model scores differently.

Answers the question a chemist actually asks — “what would I have to change to fix this?” — by searching a candidate set for the structurally closest molecule whose prediction differs by at least a given margin. Unlike an importance ranking, the result is an actionable, synthesizable alternative rather than an abstraction.

Parameters:

  • model (Any)

  • fingerprint (Any) – Used both for prediction and for the similarity search.

  • delta (float) – Minimum prediction difference for a molecule to count as a counterfactual. On a pActivity scale, 1.0 is a ten-fold change.

  • min_similarity (float) – Minimum Tanimoto similarity, so the counterfactual is a recognisable relative rather than an unrelated molecule.

  • use_proba (bool)

Examples

>>> from rdkit import Chem
>>> from sklearn.ensemble import RandomForestRegressor
>>> from qsarkit.representation import MorganFingerprint
>>> fp = MorganFingerprint(n_bits=256)
>>> library = [Chem.MolFromSmiles(s) for s in
...            ("CCO", "CCN", "CCC", "c1ccccc1", "c1ccccc1O")]
>>> model = RandomForestRegressor(n_estimators=5, random_state=0).fit(
...     fp.transform(library), [1.0, 1.1, 1.2, 8.0, 8.1]
... )
>>> explainer = CounterfactualExplainer(model, fp, delta=2.0,
...                                     min_similarity=0.0)
>>> result = explainer.explain(library[0], library)
>>> result is None or "counterfactual" in result
True

References

  • Wachter, S., Mittelstadt, B. & Russell, C. (2018). “Counterfactual Explanations without Opening the Black Box.” Harvard J. Law & Tech., 31(2), 841-887. https://doi.org/10.2139/ssrn.3063289

  • Wellawatte, G. P., Seshadri, A. & White, A. D. (2022). “Model Agnostic Generation of Counterfactual Explanations for Molecules.” Chem. Sci., 13, 3697-3705. https://doi.org/10.1039/D1SC05259D

  • Hussain, J. & Rea, C. (2010). “Computationally Efficient Algorithm to Identify Matched Molecular Pairs.” J. Chem. Inf. Model., 50(3), 339-348. https://doi.org/10.1021/ci900450m

explain(mol, candidates)[source]

Find the most similar candidate the model scores differently.

Parameters:

  • mol (Any) – The molecule to explain.

  • candidates (Sequence[Any]) – Molecules to search. Typically a virtual library or the rest of the dataset.

Returns:

counterfactual (Mol), smiles, similarity, prediction, original_prediction, delta, and transformation (the matched-pair change, when the two form a matched molecular pair). None when no candidate meets both thresholds.

Return type:

Optional[Dict[str, Any]]

References

  • Lundberg, S. M. & Lee, S.-I. (2017). “A Unified Approach to Interpreting Model Predictions.” NeurIPS 2017, 4765-4774. https://papers.nips.cc/paper/7062

  • Ribeiro, M. T., Singh, S. & Guestrin, C. (2016). “Why Should I Trust You?” KDD 2016, 1135-1144. doi:10.1145/2939672.2939778

  • Breiman, L. (2001). “Random Forests.” Machine Learning, 45(1), 5-32. doi:10.1023/A:1010933404324

  • Riniker, S. & Landrum, G. A. (2013). “Similarity Maps — A Visualization Strategy for Molecular Fingerprints and Machine-Learning Methods.” J. Cheminform., 5, 43. doi:10.1186/1758-2946-5-43

  • Polishchuk, P. et al. (2016). “Universal Approach for Structural Interpretation of QSAR/QSPR Models.” Mol. Inform., 32(9-10), 843-853. doi:10.1002/minf.201300029