Reference

• scRNA-seq: 10x H5 count files and MatrixMarket folders.
• scATAC-seq: Cell Ranger ATAC output folders with peak matrix, fragments, fragment index, and peak annotation files.
• scMultiome: Cell Ranger ARC folders containing paired RNA and ATAC output.
• Xenium: Xenium output folders with cell-feature matrix, spatial coordinates, analysis files, and optional histology alignment.
• Visium HD: Space Ranger output folders, segmented outputs, binned outputs, tissue positions, and optional histology image overlays.
• MERFISH: MERSCOPE-style cell-by-gene and cell metadata files.
• CosMx: CosMx expression matrix and metadata files.

The standard RNA workflow is provided by bakana and scran.js: QC, filtering, normalization, highly variable gene selection, PCA, UMAP, graph-based clustering, and marker detection. Cells can be filtered by detected genes, total counts, and mitochondrial percentage. Clustering uses a shared-nearest-neighbor graph and Leiden/Louvain community detection.

Entry point: imputeGeneExpression() / runSingleSamplePipeline() in src/workers/analysis.worker.js. Loads 10x HDF5 or MatrixMarket peak × cell matrices stored as SparseMatrixCSC.

Algorithm Steps

1. Drop cells with fewer than 1,500 detected peaks (QC pre-filter).
2. Find top features via empirical CDF of row sums (q5 cutoff).
3. TF-IDF normalization (Signac method 1, scaleFactor=1e4), then log1p-transformed.
4. Randomized SVD producing 50 LSI components (WASM/SIMD accelerated via randomizedSVD() in src/scatac/svd.js).
5. Drop LSI component 0 (depth-correlated); use dims 2–50.
6. Multi-sample only: Harmony batch correction on LSI embeddings.
7. UMAP on LSI or Harmony embeddings (cosine distance).
8. KNN (k=20) to SNN (Jaccard) to Louvain clustering (from src/scatac/clustering.js).

Parameters

Entry point: plotAtacGeneActivity() / plotAtacGeneActivityFromState() in src/workers/analysis.worker.js (~line 13541). Results are cached in scAtacGeneActivityCache (Map keyed by lowercase gene name).

Algorithm Steps

1. Peak lookup uses two strategies in priority order:
   1. Annotation-based (getPeaksForGene()): scan peakAnnotation rows for gene column match (case-insensitive; handles comma/semicolon multi-gene entries); resolves peak matrix row index by matching name formats chr:start-end, chr_start_end, chr-start-end.
   2. TSS fallback (getPeaksForGeneByTSS()): look up TSS in built-in reference; collect peaks whose midpoint falls within ±500 kb; chromosome names normalized (with/without chr prefix).
2. Per-cell score: sum raw peak counts across all linked peaks.
   1. Fast path: iterate CSC colPtr / rowIdx / values directly at O(nnz), with no dense allocation.
   2. Fallback: dense .column(c) access per cell.
3. Score is raw counts on a linear scale, with no log or normalization for standalone scATAC. The integration path applies log1p(count/depth × medianDepth) plus p15 background subtraction, and is triggered only when modality === 'atac-integration'.
4. Expression range clipped at p2–p98 for color scale.

Parameters

Entry point: shares the same plotAtacGeneActivity() / plotAtacGeneActivityFromState() call (~line 13541). Frontend rendering in src/components/CoveragePlot.jsx.

Algorithm Steps

1. Peak lookup via getPeaksForGene() (annotation-based); falls back to getPeaksForGeneByTSS() within 500 kb of TSS.
2. Per-cell gene activity score: sum raw peak counts across linked peaks (CSC fast path O(nnz) or dense fallback).
3. Integration normalization (atac-integration only):
   1. Depth-normalize: log1p((count / depth) × medianDepth) to correct cross-sample depth differences.
   2. Background subtraction: subtract the p15 noise floor and clip to 0, which collapses non-expressing cells and sharpens color contrast.
4. Coverage signal per cluster (Signac-style, computeCoverageByClusterRaw()):
   1. group_scale_factor = mean_depth × n_cells per cluster.
   2. signal = raw_sum / group_scale_factor × median(group_scale_factors).
   3. Fast path: single CSC scan accumulates all cluster × peak sums at once.
5. Global normalization: normalizeCoverageByGlobalMax() divides all tracks by the global max, placing all clusters on a shared [0, 1] scale.
6. Frontend rendering: clusters with the same display label are merged (weighted-average signal by cell count); values above yMax are hard-clipped; D3 area() with curveBasis (B-spline) draws smooth filled tracks. Fixed sticky header shows genomic axis (Mb), PEAKS track, and Gene track; scrollable body shows per-cluster tracks.

Parameters

Entry point: runWNNPipeline() in src/scatac/wnn.js. Combines RNA PCA and ATAC LSI modalities for multiome data. Uses buildKNN() from src/scatac/clustering.js for per-modality neighbor graphs.

Algorithm Steps

1. Build k-NN (k=20, cosine) in RNA-PCA space; z-score normalization is skipped because raw cosine correctly weights high-variance PCs.
2. Build k-NN (k=20, cosine) in ATAC-LSI space (LSI dim 0 already excluded before entry).
3. Compute per-cell modality weights from mean neighbor distance:
   w_RNA[i] = (1/avgDist_RNA) / (1/avgDist_RNA + 1/avgDist_ATAC)
4. Build proportional combined k-NN for UMAP: take round(k × w_RNA[i]) from RNA, remainder from ATAC; shared neighbors get (d_RNA + d_ATAC) / 2 (rewarded, not penalized). Calls buildWeightedCombinedKNN().
5. Run UMAP on the precomputed weighted k-NN graph with random initialization (avoids spectral-layout artifacts). Uses umap-js via setPrecomputedKNN().
6. Build a second combined k-NN for clustering, then run buildSNN() and louvain() on the combined graph.

Parameters

Entry point: runHarmony() in src/scatac/harmony.js, called from src/scatac/runMultiSamplePipeline.js (~line 141). Operates on LSI dims 2–50 (49 dimensions; component 0 dropped as depth-correlated). All internal buffers use column-major layout (d × N) for cache-efficient inner loops.

Algorithm Steps

1. L2-normalize input embeddings into Z_cos (column-major d × N, for cosine-distance clustering).
2. KMeans++ initialization (25 warmup iterations): pick first centroid randomly; select subsequent centroids proportional to squared distance from nearest existing centroid; refine with 25 standard k-means rounds.
3. Compute initial soft assignments: R[k,i] = exp(−dist[k,i] / σ[k]), normalized per cell (column-wise softmax).
4. Compute observed O[k,b] and expected E[k,b] = rowSum(R[k,:]) × Pr_b[b] batch matrices.
5. Main Harmony loop (max 20 iterations):
   1. Cluster step (max 20 k-means rounds, converge if rel. change < 1e-5 over window=3):
      1. Update centroids Y: weighted sum of Z_cos by R, L2-normalized.
      2. Distances: dist[k,i] = 2(1 − Y[:,k]·Z_cos[:,i]) (cosine via dot product).
      3. Update R with batch diversity penalty in shuffled blocks (5% cells/block): R[k,i] ∝ exp(−dist/σ) × (E/(O+E))^θ, renormalized.
   2. Correct step (ridge regression per cluster k):
      1. Solve: W = (Φ'RΦ + λI)⁻¹ Φ'RZ_orig via Gauss-Jordan inversion of (B+1)×(B+1) matrix.
      2. Zero out intercept row of W (removes global shift, keeps batch-specific correction only).
      3. Z_corr -= W × Φ × R[k,:]; L2-renormalize to get the updated Z_cos.
   3. Check Harmony convergence: objective rel. change < 1e-6.
6. Convert Z_corr from column-major back to row-major output (nCells × 49), passed directly to UMAP and clustering.

Objective function tracked per iteration:
obj = (Σ R·dist + Σ σ·R·log(R) + Σ σ·θ·O·log((O+E)/E)) × 2000/N
Terms: soft k-means error · entropy regularization · batch diversity cross-entropy penalty.

Parameters

Entry point: imputeGeneExpression() in src/workers/analysis.worker.js (~line 10951). Predicts expression of genes present in a reference scRNA-seq dataset but absent from the spatial panel (Xenium). Results are cached in imputedGeneCache (gene mapped to Float32Array) and tagged isImputed: true for visualization only; they are never used in PCA or clustering.

Algorithm Steps

1. Z-score normalize genes (column-wise) on both spatial and RNA datasets.
2. Find common genes between spatial and RNA.
3. Run PCA on each dataset (common genes only).
4. Orthogonalize PCA components via Gram-Schmidt.
5. SVD of the cross-covariance matrix to extract principal vectors.
6. Filter principal vectors by cosine similarity > 0.3.
7. Project both datasets onto RNA principal vectors.
8. Build k-NN graph (k=50, cosine distance) in projected space.
9. Predict missing genes via weighted neighbor averaging.

Input: spatial matrix from loaded Xenium data (log-normalized). scRNA reference loaded from 10x HDF5 (v2/v3) or MatrixMarket; scRNA counts normalized to 10k library size + log1p before runSpaGE().

Parameters

Entry point: runBanksyRegionSegmentation() in src/workers/analysis.worker.js (~line 11153), calling runBanksy() from src/spatial/banksy.js. Results are stored in currentResults.regionClusters, which is kept separate from transcriptomic clusters and never overwrites them.

Algorithm Steps

1. Select top 500 HVGs by per-gene variance (up to 5,000 cells sampled for speed); build dense Float32Array (nCells × nHvgs).
2. Build spatial k-NN graph (k=15) using grid-accelerated 2D search (~5–10 cells per grid cell). Reads coordinates from loadedData.spatialData; filters cells with null/NaN coords.
3. Compute Gaussian weights: w_ij = exp(−d²_ij / median(d²_i)), row-normalized per cell.
4. Compute neighbor mean matrix: N = W @ X (weighted average of spatial neighbor expression).
5. Z-score X and N independently, column-wise (per gene).
6. Assemble BANKSY matrix: [sqrt(1−λ)·zscore(X) | sqrt(λ)·zscore(N)], with shape nCells × 2·nGenes.
7. Randomized PCA on BANKSY matrix (20 components, 3 power iterations, Box-Muller normal init) with seeded RNG (seed=42).
8. Build k-NN (k=50) in BANKSY-PCA space via scran.js, running buildNeighborSearchIndex, findNearestNeighbors, buildSnnGraph, and clusterGraph in sequence. An approximate index is used if nCells > 10,000.

Parameters

Entry point: runVisiumHDSketchPipeline() in src/workers/analysis.worker.js (~line 3719), calling runSketchClustering() from src/spatial/sketchClustering.js. Triggered automatically when isVisiumHDData() returns true and nCells ≥ 50,000. Falls back to standard runClusteringAndUMAP() on any failure.

Algorithm Steps

1. Fetch PCA from bakana state (rna_pca.fetchPCs()); convert from column-major to row-major Float64Array. If fewer than 50 PCs were used, recompute via bakana.
2. Compute per-cell leverage scores: score[i] = Σ_k PC[i,k]² (L2 norm² in PCA space; rare/outlier cells score highest).
3. Two-phase sampling to select 50,000 sketch cells:
   1. Phase 1 (25% guaranteed): top 12,500 cells by leverage score taken deterministically — ensures rare populations always have representatives.
   2. Phase 2 (75% random): proportional weighted sampling from remaining pool using leverage scores as weights.
4. Fit UMAP on sketch only (k=30, minDist=0.1, Euclidean on PCA coords).
5. Cluster sketch: reuse UMAP's k-NN, truncate to snnK=20 for SNN (pruneSNN=0) and Louvain (res=1.8). Resolution is read from sketch_resolution only, not from snn_graph_cluster.
6. Split disconnected Louvain communities: BFS over the SNN subgraph per cluster, giving each contiguous component its own label.
7. Project all non-sketch cells via random-projection ANN: 8 Gaussian unit vectors in full PCA space; binary-search ±40 candidates per projection; exact PCA distances on ~640 candidates; top-5 vote determines cluster label; UMAP position = nearest sketch cell + 2.5% noise jitter.
8. Align spatial coordinates via barcode-to-coordinate mapping after clustering.

Parameters

Entry point: runTfMotifAnalysisAction() in src/workers/analysis.worker.js (~line 10512). Triggered by intent tf_motif_analysis (~line 9903); requires peakGeneLinks and genome to be present. JASPAR PWM results cached in jasparPwmCache.

Algorithm Steps

1. Find top 100 linked marker genes for the target cluster: rank peak-gene linked genes by log2FC = log2((meanIn + ε) / (meanOut + ε)) using the normalized RNA matrix (computeTopLinkedMarkers()).
2. Collect query peaks: all peaks linked to those marker genes (capped at 300).
3. Sample 150 background peaks from peaks linked to non-marker genes; supplemented from peakAnnotation if the pool is insufficient.
4. Fetch DNA sequences for all peaks from the UCSC genome browser API (https://api.genome.ucsc.edu/getData/sequence) in batches of 20 async requests. Supports hg38, hg19, mm10, mm39.
5. Load JASPAR 2024 CORE vertebrates non-redundant PWMs (~400 KB) via fetchJasparPwms() (JASPAR REST API first, bulk .txt flat file as fallback); convert PFM to PSSM with pseudocount 0.1; threshold = 60% of max possible score; discard motifs with mean IC < 0.5 bits/position.
6. Scan each peak sequence (both strands, with branch-and-bound early exit) against every PSSM via peakHitsMotif(), counting hits per motif in query vs. background.
7. Hypergeometric test per motif (hypergeomPvalue()) via log-gamma for numerical stability; Benjamini-Hochberg correction (bhCorrect()) across all motifs.
8. RENIN step 2 (for candidates with p < 0.05 and queryHits > 0):
   reninScore = max(0, log2FC) × meanExpr_cluster(TF) × Σ max(0, Pearson(TF, gene))
   summed across all marker gene expression vectors within-cluster, balancing motif specificity, TF expression level, and TF–target co-expression.
9. Deduplicate by TF name (JASPAR has multiple matrices per TF; keep highest-scoring); return top 30.

Parameters

Agent Mode uses a selected API provider to plan CellPilot tool calls. Available providers include Gemini, ChatGPT, Claude, Groq, and OpenRouter. The agent can annotate clusters, annotate selected spatial regions, plan multi-step analysis, and summarize cell-cell interaction results from CellPilot-generated tables.

• bakana: scRNA-seq analysis framework
• scran.js: WebAssembly computation (normalization, PCA, SNN, clustering)
• deck.gl: GPU-accelerated spatial visualization
• umap-js: UMAP dimensionality reduction
• Harmony: reference algorithm for multi-sample ATAC batch correction
• Seurat WNN: reference workflow for weighted nearest neighbor multimodal integration
• BANKSY: reference algorithm for spatial domain segmentation
• SpaGE: reference algorithm for spatial gene expression imputation
• JASPAR 2024: transcription factor binding motif database

CellPilot provides a self-contained MCP (Model Context Protocol) server that exposes its analysis tools to MCP-compatible desktop assistants, currently Claude Desktop and Codex. The server bundles the analysis engine, so the separate CellPilot desktop app is not required. The assistant calls the server's tools to load data, run workflows, and return plots and result tables.

• Distribution: macOS bundles (Apple Silicon and Intel), .zip wrapped. Claude uses a .mcpb bundle; Codex uses a .dmg installer. Both are double-click installs.
• Connection: the assistant communicates with CellPilot over MCP; the assistant must be the installed desktop application, not the web version.
• Relationship to Agent Mode: Agent Mode runs inside CellPilot with a user-supplied API key, while the MCP server lets an external assistant call CellPilot's tools.

Setup steps are in the Connect via MCP tutorial; downloads are on the Download page.