scAdam model explainability

scAdam model explainability#

Model explainability refers to the ability to understand and articulate how a machine learning model makes its predictions or decisions.

scNoah uses SHAP (SHapley Additive exPlanations) for model explainability analysis.

SHAP values are derived from cooperative game theory, specifically using Shapley values, which quantify the contribution of each player (or feature/gene) to the total payoff (or prediction). This approach allows for a fair distribution of contributions among features, similar to how players in a coalition share the rewards based on their input.

The SHAP framework operates under the principle that the prediction for an instance can be expressed as the sum of the average prediction and the contributions of each feature. Mathematically, this can be represented as:

\[prediction = f_0 + \sum_{i=1}^M f_i\]

where \(prediction\) is the model’s prediction, \(f_0\) is the expected value of the output, \(f_i\) are the SHAP values for each feature.

SHAP provides a framework for both local explainability and global explainability in machine learning models.

  1. Local Explainability: Focuses on explaining individual predictions (cells) made by the model. It seeks to clarify why a specific output was generated for a given input.

  2. Global Explainability: Aims to provide an overall understanding of the model’s behavior across all predictions (cells), identifying the main factors that influence outputs.

[1]:

# Python packages
import warnings
warnings.simplefilter('ignore')

import scanpy as sc
import scparadise
import shap
import os
import matplotlib.pyplot as plt
import numpy as np
import pandas as pd

sc.set_figure_params(dpi = 120)

Recommendations about dataset#

Our models trained using shifted logarithm normalized data. We recommend shifted logarithm data normalization method to proper usage of our models: sc.pp.normalize_total(adata, target_sum=None) sc.pp.log1p(adata) adata.raw = adata

[2]:

# Load dataset from 10x Genomics
url = "https://cf.10xgenomics.com/samples/cell-exp/6.1.0/10k_PBMC_3p_nextgem_Chromium_X/10k_PBMC_3p_nextgem_Chromium_X_filtered_feature_bc_matrix.h5"
adata = sc.read_10x_h5("dataset.h5", backup_url = url)
adata.var_names_make_unique()
adata

100%|█████████████████████████████████████| 37.7M/37.7M [00:00<00:00, 78.4MB/s]

[2]:

AnnData object with n_obs × n_vars = 11996 × 36601
    var: 'gene_ids', 'feature_types', 'genome'

QC#

[3]:

# mitochondrial genes, "MT-" for human, "Mt-" for mouse
adata.var["mt"] = adata.var_names.str.startswith("MT-")
# ribosomal genes
adata.var["ribo"] = adata.var_names.str.startswith(("RPS", "RPL"))
# hemoglobin genes
adata.var["hb"] = adata.var_names.str.contains("^HB[^(P)]")

[4]:

sc.pp.calculate_qc_metrics(
    adata, qc_vars=["mt", "ribo", "hb"], inplace=True, log1p=True
)

[5]:

sc.pl.violin(
    adata,
    ["n_genes_by_counts", "total_counts", "pct_counts_mt", 'pct_counts_ribo', 'pct_counts_hb'],
    jitter=0.4,
    multi_panel=True,
)
../../../_images/tutorials_notebooks_scNoah_scAdam_model_explainability_8_0.png 
[6]:

sc.pp.filter_cells(adata, min_genes=200)
sc.pp.filter_genes(adata, min_cells=3)

[7]:

# Detect doublets
sc.pp.scrublet(adata)

[8]:

# Remove doublets + other QC metrics
adata = adata[adata.obs['predicted_doublet'] == False]
sc.pp.filter_cells(adata, max_genes = 5000)
sc.pp.filter_cells(adata, max_counts = 20000)
adata = adata[adata.obs['pct_counts_mt'] < 15]

Normalization, HVG, neighbors, PCA, UMAP#

We recommend using shifted logarithm data normalization as described here.

[9]:

# Saving count data
adata.layers["counts"] = adata.X.copy()
# Normalization (shifted logarithm)
sc.pp.normalize_total(adata, target_sum=None)
sc.pp.log1p(adata)
# scParadise use normalized data in adata.raw!!!
adata.raw = adata
# HVG
sc.pp.highly_variable_genes(adata, n_top_genes=2000)
# PCA
sc.tl.pca(adata)
# Nearest neighbors analysis
sc.pp.neighbors(adata, n_neighbors=10, n_pcs=20)
# UMAP
sc.tl.umap(adata)

Cell type prediction using scAdam#

[10]:

# Available models for cell type annotation
df = scparadise.scadam.available_models(query = 'PBMC')
df

[10]:
Model name Tissue/Organ Species Description Suspension Accuracy Balanced Accuracy Number of Levels
0 Human_PBMC blood human Peripheral blood mononuclear cells of healthy ... cells 0.979 0.979 3 
[11]:

# Download model for cell type prediction
scparadise.scadam.download_model('Human_PBMC', save_path='')

[12]:

# Predict cell types using trained model
adata = scparadise.scadam.predict(
    adata,
    path_model = 'Human_PBMC_scAdam'
)

Cell type annotation using the 2nd version of the scAdam model
scAdam model with unknown detector loaded from Human_PBMC_scAdam
Gene alignment:
  Model features: 947
  Matched features: 935 (98.7%)

Predicting: 100%|█████████████████████████████| 43/43 [00:00<00:00, 103.71it/s]

Added cell type column: pred_celltype_l1
Added probabilities column: pred_celltype_l1_probability
Added cell type column: pred_celltype_l2
Added probabilities column: pred_celltype_l2_probability
Added cell type column: pred_celltype_l3
Added probabilities column: pred_celltype_l3_probability

[14]:

# Visualise predicted cell types levels and prediction probabilities
sc.set_figure_params(figsize=[8, 5])
sc.pl.umap(
    adata,
    color = [
        'pred_celltype_l1',
        'pred_celltype_l1_probability',
        'pred_celltype_l2',
        'pred_celltype_l2_probability',
        'pred_celltype_l3',
        'pred_celltype_l3_probability'
    ],
    frameon = False,
    legend_loc = 'on data',
    legend_fontsize = 12,
    legend_fontoutline = 2,
    ncols = 2,
    size = 30
)
../../../_images/tutorials_notebooks_scNoah_scAdam_model_explainability_18_0.png

Create explanations matrix#

Compute SHAP values for a specific cell type predicted using scAdam model.

[15]:

# Compute explanations for CD8+ T cells
random_state = 0
num_cells = 100
explanations  = scparadise.scnoah.explain(
    adata,
    celltype = 'CD8 T',
    path_model = 'Human_PBMC_scAdam',
    num_cells = num_cells,
    random_state = random_state
)

scAdam model with unknown detector loaded from Human_PBMC_scAdam

PermutationExplainer explainer: 101it [01:59,  1.25s/it]

SHAP explanations for cell type 'CD8 T' were completed for 100 cells.



[20]:

# General explainability
# Calculate gene impoortances
gene_importances = scparadise.scnoah.feature_importance_df(
    explanations
)
# Show 15 most important genes for CD8+ T cells annotation
gene_importances.head(12)

[20]:
gene_name gene_importance mean_shap_value
0 CCL5 0.045810 0.002262
1 CD8B 0.034693 -0.000737
2 TXK 0.032035 -0.000124
3 GZMA 0.031599 0.000551
4 CD8A 0.028938 -0.001154
5 LINC02446 0.026669 -0.001247
6 ALOX5AP 0.025044 -0.000834
7 CORO1B 0.024495 0.001444
8 NUCB2 0.022090 0.001117
9 AQP3 0.021154 0.000332
10 CLEC2D 0.019695 0.000417
11 LINC00861 0.018677 -0.000530 
[25]:

# Visualization of gene importancies on cell embeddings
sc.pl.umap(
    adata,
    color = gene_importances.head(11)['gene_name'].tolist() + ['pred_celltype_l2'],
    frameon = False,
    cmap = 'coolwarm',
    vmax = 'p99',
    size = 30,
    ncols = 4
)
../../../_images/tutorials_notebooks_scNoah_scAdam_model_explainability_22_0.png

gene_importances highlights the importance of CD8A and CD8B for annotating CD8+ T cells, which is biologically expected because both genes are canonical CD8+ T-cell markers.

Several highly ranked genes are also associated with distinct CD8+ T-cell states. GZMK, GZMA, CCL5, NKG7, and KLRD1 are commonly expressed by memory and cytotoxic/effector CD8+ T-cell populations, whereas IL7R and CCR7 are characteristic of naïve and central-memory-like CD8+ T cells. Thus, the gene-importance results capture both the core CD8+ T-cell identity and transcriptional differences among CD8+ T-cell subtypes.

Analysis of explanations#

  1. Positive SHAP value indicates that a feature (gene) increases the model’s output for the cell type being explained, relative to the model’s baseline prediction.

  2. Negative SHAP value indicates that the feature (gene) decreases that output.

In other words, for a specific cell and target cell type, a gene with a positive SHAP value supports the predicted label, whereas a gene with a negative SHAP value argues against it. The magnitude of the SHAP value reflects the strength of that contribution—not necessarily the gene’s overall biological importance. SHAP values are local, instance-specific contributions that collectively explain the difference between an individual prediction and the model’s expected baseline output.

[29]:

# Visualize feature explanations.
# The bar plot shows the genes with the greatest overall importance
# for CD8+ T-cell annotation.
plt.grid(visible=False, axis="both")
shap.plots.bar(explanations, max_display=12)
../../../_images/tutorials_notebooks_scNoah_scAdam_model_explainability_25_0.png 
[30]:

# Visualize SHAP values for individual cells.
shap.plots.heatmap(explanations, max_display=15, show=False)
plt.grid(visible=False, axis="both")
../../../_images/tutorials_notebooks_scNoah_scAdam_model_explainability_26_0.png

The bar plot shows the genes with the greatest overall contribution to CD8+ T-cell annotation. The heatmap displays SHAP values for individual cells, allowing inspection of how each gene contributes to the predicted CD8+ T-cell label for each cell.

Difference between CD8+ T cells#

Differences in SHAP values suggest the presence of at least two distinct groups of CD8+ T cells:

  1. Cells with positive SHAP values for CCL5 and negative SHAP values for GZMK.

  2. Cells with positive or zero SHAP values for GZMK and negative or zero SHAP values for CCL5.

For individual cells, SHAP values can be displayed together with gene-expression levels.

Warning! shap.plots.heatmap() performs hierarchical clustering automatically. To identify the cell corresponding to a position in the heatmap, perform the same clustering on explanations: cluster_order = explanations.hclust() You can then use the instance number displayed by shap.plots.heatmap() together with cluster_order to identify the corresponding cell in the original explanations object.

[45]:

# Memory CD8+ T cell
# The numbers shown to the left of the gene names represent gene-expression values.
# `explanations.hclust()[1]` maps the second cell in the heatmap ordering
# to its index in the original `explanations` object.
cell_index = explanations.hclust()[50]

shap.plots.bar(
    explanations[cell_index],
    max_display=15,
    show=False,
)
plt.grid(visible=False, axis="both")

# Index of the corresponding cell in the original explanations object
cell_index

[45]:

np.int32(29)
../../../_images/tutorials_notebooks_scNoah_scAdam_model_explainability_30_1.png 
[38]:

# Naïve CD8+ T cell
# The numbers shown to the left of the gene names represent gene-expression values.
# `explanations.hclust()[90]` maps the 91st cell in the heatmap ordering
# to its index in the original `explanations` object.
cell_index = explanations.hclust()[90]

shap.plots.bar(
    explanations[cell_index],
    max_display=15,
    show=False,
)
plt.grid(visible=False, axis="both")

# Index of the corresponding cell in the original explanations object
cell_index

[38]:

np.int32(62)
../../../_images/tutorials_notebooks_scNoah_scAdam_model_explainability_31_1.png

The left-hand numbers next to gene names represent the gene-expression values for the selected cell. The bars show the SHAP values: positive values increase the model output for the explained CD8+ T-cell prediction, whereas negative values decrease it.

In this example, explanations.hclust()[1] retrieves the index of the second cell in the hierarchical-clustering order used by the heatmap. This index identifies the corresponding cell in the original explanations object.

[39]:

# Extract the same 100 cells from a dataset
adata_subset = adata[adata.obs['pred_celltype_l2'] == 'CD8 T'].copy()
indexes = pd.DataFrame(index = adata_subset.obs_names).sample(n = num_cells, axis = 0, random_state = random_state).index.tolist()
adata_subset = adata_subset[indexes]

[50]:

# Visualization of gene importances on cell embeddings
# CD8 TCM - central memory
# CD8 TEM - effector memory
sc.pl.embedding(
    adata_subset,
    color = [
        'CHST12',
        'GNLY',
        'CCL5',
        'pred_celltype_l3'
    ],
    basis = 'X_umap',
    frameon = False,
    legend_loc = 'right margin',
    palette = 'tab20',
    vmax = 'p99',
    cmap = 'coolwarm',
    ncols = 3,
    size = 500
)
../../../_images/tutorials_notebooks_scNoah_scAdam_model_explainability_34_0.png

Marker gene search based on gene importances#

[51]:

# Compute explanations for plasmacytoid dendritic cells (pDC)
explanations  = scparadise.scnoah.explain(
    adata,
    celltype = 'pDC',
    num_cells = 100,
    path_model = 'Human_PBMC_scAdam'
)

scAdam model with unknown detector loaded from Human_PBMC_scAdam

PermutationExplainer explainer: 88it [02:22,  1.72s/it]

SHAP explanations for cell type 'pDC' were completed for 87 cells.



[52]:

# Calculate gene impoortances
gene_importances = scparadise.scnoah.feature_importance_df(
    explanations
)
# Show 10 most important genes for pDC cell annotation
gene_importances.head(12)

[52]:
gene_name gene_importance mean_shap_value
0 NUCB2 0.005982 0.000023
1 OPN3 0.005191 0.000031
2 CCDC88A 0.004623 -0.000168
3 PTGDS 0.004364 -0.000183
4 PACSIN1 0.003192 -0.000103
5 CD68 0.003032 -0.000112
6 TXNDC5 0.002829 0.000103
7 LILRB4 0.002777 0.000037
8 BEX4 0.002571 -0.000037
9 CLIC3 0.002479 0.000107
10 MYO1E 0.002378 0.000028
11 AFF3 0.002320 0.000078 
[53]:

# Visualization of gene importancies on cell embeddings
sc.pl.umap(
    adata,
    color = gene_importances.head(11)['gene_name'].tolist() + ['pred_celltype_l2'],
    frameon = False,
    cmap = 'coolwarm',
    vmax = 'p99',
    size = 30,
    ncols = 4
)
../../../_images/tutorials_notebooks_scNoah_scAdam_model_explainability_38_0.png

PACSIN1 and PTGDS support a pDC-like identity, while LILRB4 provides additional, although less specific, support.

[54]:

# Save anndata with predicted annotations
adata.write_h5ad('adata_predicted.h5ad')

[55]:

import session_info2
session_info2.session_info()

[55]:
PackageVersion
pandas2.3.3
numpy2.4.6
scanpy1.11.5
scparadise1.2.0
shap0.49.1
matplotlib3.11.0
anndata0.12.16
ComponentInfo
Python3.11.15 (main, Jun 11 2026, 15:20:16) [GCC 14.3.0]
OSLinux-6.8.0-138-generic-x86_64-with-glibc2.39
Updated2026-09-04 10:32
Dependencies
DependencyVersion
narwhals2.22.1
pynndescent0.6.0
scverse-misc0.0.3
opt_einsum3.4.0
pillow12.2.0
wcwidth0.8.1
rich15.0.0
fonttools4.63.0
session-info20.4.1
texttable1.7.0
statsmodels0.14.6
igraph1.0.0
decorator5.3.1
ipykernel7.3.0
cloudpickle3.1.2
threadpoolctl3.6.0
kiwisolver1.5.0
aiohttp3.14.1
imbalanced-learn0.14.2
sparse0.18.0
platformdirs4.10.0
typing_extensions4.15.0
cycler0.12.1
prompt_toolkit3.0.52
argparse1.4.0 (1.1)
joblib1.5.3
dataclasses0.6
h5py3.16.0
aiosignal1.4.0
pyarrow24.0.0
legacy-api-wrap1.5
comm0.2.3
PyYAML6.0.3
cuda-bindings13.3.1
scikit-image0.26.0
setuptools81.0.0
numba0.65.1
aiohappyeyeballs2.6.2
pyparsing3.3.2
requests2.34.2
parso0.8.7
scikit-learn1.9.0
packaging26.0
fsspec2026.4.0
google-crc32c1.8.0
pytz2026.2
Cython3.3.0
scipy1.17.1
executing2.2.1
mudata0.3.8
urllib32.7.0
pytorch-tabnet4.1.0
certifi2026.5.20 (2026.05.20)
traitlets5.15.1
jedi0.20.0
pathlib1.0.1
tornado6.5.7
lazy-loader0.5
seaborn0.13.2
propcache0.5.2
asttokens3.0.1
yarl1.24.2
stack-data0.6.3
umap-learn0.5.12
defusedxml0.7.1
leidenalg0.12.0
llvmlite0.47.0
charset-normalizer3.4.7
donfig0.8.1.post1
python-dateutil2.9.0.post0
xarray2026.4.0
six1.17.0
pycparser3.0 (3.00)
jupyter_core5.9.1
natsort8.4.0
torch2.12.0 (2.12.0+cu130)
Pygments2.20.0
sklearn-compat0.1.6
plottable0.1.5
muon0.1.7
slicer0.0.8
attrs26.1.0
optuna4.9.0
pyzmq27.1.0
numcodecs0.16.5
psutil7.2.2
colorlog6.10.1
idna3.18
matplotlib-inline0.2.2
patsy1.0.2
debugpy1.8.21
cuda-pathfinder1.5.5
multidict6.7.1
jupyter_client8.9.1
cffi2.0.0
tqdm4.68.2
ipython9.14.1
zarr3.1.6
frozenlist1.8.0
typing3.7.4.3
pure_eval0.2.3
Copyable Markdown 
| Package    | Version |
| ---------- | ------- |
| pandas     | 2.3.3   |
| numpy      | 2.4.6   |
| scanpy     | 1.11.5  |
| scparadise | 1.2.0   |
| shap       | 0.49.1  |
| matplotlib | 3.11.0  |
| anndata    | 0.12.16 |

| Dependency         | Version                |
| ------------------ | ---------------------- |
| narwhals           | 2.22.1                 |
| pynndescent        | 0.6.0                  |
| scverse-misc       | 0.0.3                  |
| opt_einsum         | 3.4.0                  |
| pillow             | 12.2.0                 |
| wcwidth            | 0.8.1                  |
| rich               | 15.0.0                 |
| fonttools          | 4.63.0                 |
| session-info2      | 0.4.1                  |
| texttable          | 1.7.0                  |
| statsmodels        | 0.14.6                 |
| igraph             | 1.0.0                  |
| decorator          | 5.3.1                  |
| ipykernel          | 7.3.0                  |
| cloudpickle        | 3.1.2                  |
| threadpoolctl      | 3.6.0                  |
| kiwisolver         | 1.5.0                  |
| aiohttp            | 3.14.1                 |
| imbalanced-learn   | 0.14.2                 |
| sparse             | 0.18.0                 |
| platformdirs       | 4.10.0                 |
| typing_extensions  | 4.15.0                 |
| cycler             | 0.12.1                 |
| prompt_toolkit     | 3.0.52                 |
| argparse           | 1.4.0 (1.1)            |
| joblib             | 1.5.3                  |
| dataclasses        | 0.6                    |
| h5py               | 3.16.0                 |
| aiosignal          | 1.4.0                  |
| pyarrow            | 24.0.0                 |
| legacy-api-wrap    | 1.5                    |
| comm               | 0.2.3                  |
| PyYAML             | 6.0.3                  |
| cuda-bindings      | 13.3.1                 |
| scikit-image       | 0.26.0                 |
| setuptools         | 81.0.0                 |
| numba              | 0.65.1                 |
| aiohappyeyeballs   | 2.6.2                  |
| pyparsing          | 3.3.2                  |
| requests           | 2.34.2                 |
| parso              | 0.8.7                  |
| scikit-learn       | 1.9.0                  |
| packaging          | 26.0                   |
| fsspec             | 2026.4.0               |
| google-crc32c      | 1.8.0                  |
| pytz               | 2026.2                 |
| Cython             | 3.3.0                  |
| scipy              | 1.17.1                 |
| executing          | 2.2.1                  |
| mudata             | 0.3.8                  |
| urllib3            | 2.7.0                  |
| pytorch-tabnet     | 4.1.0                  |
| certifi            | 2026.5.20 (2026.05.20) |
| traitlets          | 5.15.1                 |
| jedi               | 0.20.0                 |
| pathlib            | 1.0.1                  |
| tornado            | 6.5.7                  |
| lazy-loader        | 0.5                    |
| seaborn            | 0.13.2                 |
| propcache          | 0.5.2                  |
| asttokens          | 3.0.1                  |
| yarl               | 1.24.2                 |
| stack-data         | 0.6.3                  |
| umap-learn         | 0.5.12                 |
| defusedxml         | 0.7.1                  |
| leidenalg          | 0.12.0                 |
| llvmlite           | 0.47.0                 |
| charset-normalizer | 3.4.7                  |
| donfig             | 0.8.1.post1            |
| python-dateutil    | 2.9.0.post0            |
| xarray             | 2026.4.0               |
| six                | 1.17.0                 |
| pycparser          | 3.0 (3.00)             |
| jupyter_core       | 5.9.1                  |
| natsort            | 8.4.0                  |
| torch              | 2.12.0 (2.12.0+cu130)  |
| Pygments           | 2.20.0                 |
| sklearn-compat     | 0.1.6                  |
| plottable          | 0.1.5                  |
| muon               | 0.1.7                  |
| slicer             | 0.0.8                  |
| attrs              | 26.1.0                 |
| optuna             | 4.9.0                  |
| pyzmq              | 27.1.0                 |
| numcodecs          | 0.16.5                 |
| psutil             | 7.2.2                  |
| colorlog           | 6.10.1                 |
| idna               | 3.18                   |
| matplotlib-inline  | 0.2.2                  |
| patsy              | 1.0.2                  |
| debugpy            | 1.8.21                 |
| cuda-pathfinder    | 1.5.5                  |
| multidict          | 6.7.1                  |
| jupyter_client     | 8.9.1                  |
| cffi               | 2.0.0                  |
| tqdm               | 4.68.2                 |
| ipython            | 9.14.1                 |
| zarr               | 3.1.6                  |
| frozenlist         | 1.8.0                  |
| typing             | 3.7.4.3                |
| pure_eval          | 0.2.3                  |

| Component | Info                                               |
| --------- | -------------------------------------------------- |
| Python    | 3.11.15 (main, Jun 11 2026, 15:20:16) [GCC 14.3.0] |
| OS        | Linux-6.8.0-138-generic-x86_64-with-glibc2.39      |
| Updated   | 2026-09-04 10:32                                   |