2.1 Load the DepMap reference data

load_depmap() downloads and caches the DepMap bulk expression, single-cell expression, drug response (AUC), and cell line annotations. The first run downloads about 500 MB; afterwards the cache is reused.

load_depmap(read = TRUE, mirror = TRUE)

The mirror = TRUE argument uses GitHub mirrors for faster downloads in regions with limited connectivity. You can also add custom mirrors:

add_mirrors("https://my-mirror.example.com/PERCEPTIONx")

2.2 Load pre-trained models

PERCEPTIONx ships 44 pre-trained models for FDA-approved drugs, so you can predict without training:

# Load a single model
models <- load_model("abemaciclib", read = TRUE)

# Or several at once
models <- load_model(c("abemaciclib", "erlotinib"), read = TRUE)

If you want to train models for other drugs or settings instead, see section 3.1.

2.3 Prepare the clone-level expression data

PERCEPTIONx models are trained on rank-normalized expression, and clone-level prediction needs clones defined first. prepare_data() does both: it clusters the cells with Seurat (UMAP or t-SNE), merges each clone’s expression, and rank-normalizes it — one call:

prepared <- prepare_data(my_expression_matrix, method = "umap")

The returned list carries everything downstream needs: clone_expression_rnorm (the rank-normalized clone-level matrix) and clone_counts (clone abundances per patient). Rank normalization happens inside prepare_data() — it calls rank_normalization_mat() on the merged clone expression — so you do not run it yourself on this path. See section 3.5 for the details and the skip_clustering option.

2.4 Predict drug sensitivity

Prediction runs in two stages: clone-level scoring, then patient-level aggregation.

clone_pred <- predict_drugs(models, prepared$clone_expression_rnorm)

patient_pred <- predict_patients(clone_pred, prepared)

2.5 Visualize

Each plotting function takes the corresponding prediction output and returns a ggplot. This example overlays predicted viability on a t-SNE embedding:

plot_tsne_response(
  tsne_data = tsne_data,
  color_var = "viability_scaled",
  title     = "Drug Response"
)

The full set of plot functions — clone distribution, lollipop, ROC, boxplot, UMAP overlays, and a composite panel — is listed in section 3.4.

3.1 Train your own models

Training is only needed when you want drugs or settings not covered by the 44 pre-trained models. train_models() presets every argument — even the drug list defaults to all 44 supported drugs — so the whole call can be one line:

models <- train_models(drug_list = "erlotinib")

That is all it takes — cancer type (PanCan), model type (glmnet), cross-validation, and feature ranking are handled with sensible defaults. The knobs you might actually touch:

One thing worth knowing: output_dir defaults to "./models", meaning every trained model is written to disk as a single timestamped RDS file there (and also returned in the models object).

The returned models object plugs straight into predict_drugs() (section 2.4), and can be saved or re-loaded later:

saveRDS(models, "my_models.RDS")
models <- readRDS("my_models.RDS")

3.2 Evaluate model performance

Two helpers answer how well the trained models stratify responders from non-responders:

perf <- compare_performance(models, threshold = 0.3, verbose = TRUE)
plot_model_performance(models)

sig_models <- get_significant_models(
  models,
  min_correlation = 0.3,
  max_pvalue = 0.05
)

# One patient index at a time; returns the patient's pseudo-bulk vector.
# The clone-level expression matrix and the abundance table come straight
# out of prepare_data() (clone columns named "Patient@@Clone").
pseudo_bulk_list <- lapply(
  seq_len(nrow(prepared$clone_counts)),
  function(x) each_patient_pseudo_bulk(
    x = x,
    comb_viability_df = NULL,               # retained for compatibility, unused
    Clone_Counts_per_patients = prepared$clone_counts,
    clone_Level_z_expression_df = prepared$clone_expression_rnorm
  )
)

3.3 Validate predictions on DepMap cell lines

This is a verification exercise, not part of the prediction path: it predicts viability for held-out DepMap cell lines and checks the scores against known drug response. It is the pattern used to produce the validation figures in the paper.

cellLines_test <- get_cellLine_list(
  infunc_cancerType = "PanCan",
  infunc_drugName   = "erlotinib",
  exclude_cancer    = "PanCan",
  infunc_response   = get_response_matrix("erlotinib")
)[[2]]

test_cells <- DepMap$metadata_CPM_scRNA$NAME[
  DepMap$metadata_CPM_scRNA$DepMap_ID %in% cellLines_test]
expr_test <- DepMap$CPM_scRNA_CCLE_rnorm[, test_cells, drop = FALSE]

cell_pred <- predict_drugs(models, expr_test)

For your own patient data, skip this block entirely and go straight to predict_drugs() as in section 2.4.

3.4 The plotting suite

Several plotting functions accept tooltip = TRUE (default): when the ggiraph package is installed, points and bars get hover tooltips (clone id, viability score, proportion, FPR/TPR). This applies to plot_tsne_response(), plot_clone_umap(), plot_clone_distribution(), plot_clone_viability(), plot_roc_curve(), plot_response_boxplot(), and plot_model_performance() — the remaining functions (plot_model_roc(), plot_patient_response_panel(), plot_tsne_biomarker_viability(), plot_seurat_clustering()) do not take a tooltip argument. Set tooltip = FALSE for a plain static ggplot with the identical layout. Tooltip text is auto-built from the plot data; plot_tsne_response(), plot_clone_umap(), plot_clone_distribution(), and plot_clone_viability() also accept a custom tooltip_col column.

# t-SNE / UMAP with drug response overlay
plot_tsne_response(tsne_data, color_var = "viability_scaled", title = "Drug Response")

# Biomarker vs. viability side-by-side
plot_tsne_biomarker_viability(
  tsne_data, biomarker_var = "biomarker_scaled", viability_var = "viability_scaled",
  biomarker_label = "Biomarker Expression", viability_label = "Drug Viability"
)

# Clone abundance stacked bar chart
plot_clone_distribution(clone_distribution = clone_dist_df, response_var = "response")

# Clone-level viability lollipop
plot_clone_viability(
  clone_viability = clone_viability_df,
  viability_var = "comb_viability",
  weights_var = "weights"          # optional: point size by clone proportion
)

# Responder vs. non-responder boxplot
plot_response_boxplot(
  exp_vs_pred = exp_vs_pred,
  response_var = "response",
  predicted_var = "predicted_viability"
)

# ROC curve with AUC
plot_roc_curve(response = response, predictor = predictor, smooth_curve = TRUE)

# Seurat clustering and UMAP
result <- plot_seurat_clustering(expr_matrix)
print(result$umap_plot)

# Composite patient response panel
plot_patient_response_panel(
  clone_distribution = clone_dist_df,
  clone_viability = clone_viability_df,
  exp_vs_pred = exp_vs_pred,
  viability_col = "comb_viability"
)

Interactive example:

plot_clone_viability(
  clone_viability = clone_viability_df,
  viability_var = "comb_viability",
  weights_var = "weights",
  tooltip = TRUE
)

3.5 prepare_data() and Seurat clustering

The prepare_data() function turns a raw expression matrix into the clone-level inputs used by prediction and visualization: it clusters cells with Seurat (UMAP or t-SNE), assigns clones, merges each clone’s expression, and rank-normalizes it (internally calling rank_normalization_mat()), returning a list with clone_expression_rnorm and clone_counts.

prepared <- prepare_data(expr_matrix, method = "umap", seurat_resolution = 0.5)
clone_pred <- predict_drugs(models, prepared$clone_expression_rnorm)
patient_pred <- predict_patients(clone_pred, prepared)

Useful options: seurat_resolution (Seurat clustering resolution), genes_to_use (restrict the gene space before clustering), and skip_clustering = TRUE — if your data is already clone-level (one column per clone), prepare_data() skips clustering and only builds the normalized clone matrix and counts.

Clone identity is carried in the row names as Patient@@Clone; parse_clone_keys() splits them, and clone_mean_expression() / build_clone_counts() work on the cloned data for downstream analyses.

3.6 The Shiny web application

PERCEPTION-shiny wraps the whole pipeline in a point-and-click dashboard. Heavy computation (training, Seurat clustering, prediction, plotting) runs in background worker processes, so the interface stays responsive even while a large job runs. See the Shiny App Guide for a full walkthrough.

library(PERCEPTIONx)
run_perception_app()          # starts the app in your browser

or directly from the source tree:

shiny::runApp(system.file("shiny", "app", package = "PERCEPTIONx"))

The app has six tabs — Home, Data, Train, Predict, Visualize, and Help — and ships with a Load Demo button that generates a small synthetic dataset (49 genes x 400 cells x 20 patients) to smoke-test the whole flow. Figures in the app are interactive SVG (ggiraph-based), with static downloads available as 600 dpi PNG or vector PDF/SVG.

3.7 End-to-end pipeline script

A ready-to-run script that executes the full workflow (train -> predict -> plot) and saves all output figures lives in the source repository under tools/ (not shipped with the installed package). After cloning the repo:

source("tools/test_pipeline.R")

3.8 Agent skills

The repository ships two agent skills (.claude/skills/) that let AI coding assistants drive the package correctly instead of guessing its API:

• perceptionx-pipeline — the full workflow from this tutorial: load reference data, prepare patient data, predict, and evaluate.
• perceptionx-result — interpretation of prediction tables, model performance, and the plotting suite.

Claude Code picks them up automatically. For other agents (e.g. Codex), copy them into that agent’s skills directory:

cp -r .claude/skills/perceptionx-pipeline ~/.codex/skills/
cp -r .claude/skills/perceptionx-result ~/.codex/skills/