piaso.tl — tools

Normalization, dimensionality reduction, clustering and annotation.

Function	What it does
`analyzeMarkers`	Analyze gene list(s) to infer potential cell types.
`calculateScoreParallel`	Compute gene set scores in parallel using shared memory for efficiency.
`calculateScoreParallel_multiBatch`	Calculate gene set scores for each adata batch in parallel using shared memory. Different marker gene sets will be calculated in parallel as well.
`compute_tfidf_stats`	Compute TF-IDF statistics without materializing the full matrix.
`cospecificity_trans`	Trans co-specificity between candidate feature pairs (e.g. TF→target).
`ensure_infog_params`	Return cached `{modality}_infog_params` (or compute now and cache). Falls through to legacy `infog_params` once with a DeprecationWarning when only the unprefixed key is present. A cached payload whose `cell_depth` length no longer matches `ds.n_cells` is treated as a miss (self-heals a cytome filtered before the cached-stats invalidation fix).
`ensure_tfidf_params`	Return cached `{modality}_tfidf_params` (or compute now and cache).
`getMarkers`	Query PIASOmarkerDB for cell type marker genes.
`infog`	INFOG normalization of single-cell RNA sequencing data.
`infog_svd`	INFOG normalization → HVG selection → SVD in one call.
`leiden`	Leiden clustering of the neighbour graph.
`leiden_local`	Leiden clustering within each group of a coarse clustering.
`neighbors`	Build kNN graph and compute fuzzy simplicial set connectivities.
`PIASOmarkerDB`	Python client for accessing PIASOmarkerDB.
`predictCellTypeByGDR`	Predicts cell types in a query dataset (`adata`) using the GDR dimensionality reduction method based on a reference dataset (`adata_ref`). To use GDR for dimensionality reduction, please refer to `piaso.tl.runGDR` or `piaso.tl.runGDRParallel`.
`predictCellTypeByMarker`	Predict cell types using marker genes and optionally smooth predictions.
`projectGDR`	Project `data` into the frozen GDR space of `reference`.
`queryPIASOmarkerDB`	Query PIASOmarkerDB for cell type marker genes.
`read_selected_mask`	The boolean ‘highly_variable’ column of a modality’s feature table.
`run_TFIDF`	Compute TF-IDF normalization for peak count data.
`runCOSGParallel`	Run COSG on batches in parallel using shared memory and multiprocessing.
`runGDR`	Run GDR (marker Gene-guided dimensionality reduction) on single-cell data.
`runGDRParallel`	Run GDR (marker Gene-guided dimensionality reduction) on single-cell data.
`runHarmony`	Run Harmony batch correction on an embedding.
`runSCALAR`	Calculates ligand-receptor interaction scores, computes permutation-based p-values using a vectorized approach, and corrects for multiple testing using FDR for each cell type-cell type pair independently.
`runSVD`	Truncated SVD dimensionality reduction.
`runSVDLazy`	Deprecated alias for :func:`infog_svd`. Use `piaso.tl.infog_svd()` instead.
`score`	Compute gene-set enrichment scores for one or more gene sets, on AnnData or a Cytome dataset.
`smoothCellTypePrediction`	Smooth cell type predictions using k-nearest neighbors in a low-dimensional embedding.
`specificity_matrix`	COSG λ specificity matrix `S` [features × cell_types], for AnnData or cytome. Cytome reuses the cached :func:`_ensure_cosg_dense`; AnnData runs `cosg.cosg` at full `n_genes_user` and pivots `uns` into a dense frame.
`stitchSpace`	Performs a batch correction using a BBKNN graph that has been pruned based on marker gene overlap between batch-specific clusters. Overlap check uses local markers and optionally global markers (controlled by `filter_use_global_markers`).
`umap`	Compute UMAP embedding from precomputed kNN graph.

analyzeMarkers

analyzeMarkers(
    genes: "Union[List[str], 'pd.DataFrame', Dict[str, List[str]]]",
    n_top_genes: 'int' = 50,
    species: 'str' = None,
    tissue: 'str' = None,
    studies: 'Union[str, List[str]]' = None,
    min_genes: 'int' = 1,
    exclude_cell_types: 'List[str]' = None,
    exclude_studies: 'List[str]' = None,
)

Signature defaults

genes, n_top_genes=50, species=None, tissue=None, studies=None, min_genes=1, exclude_cell_types=None, exclude_studies=None

Analyze gene list(s) to infer potential cell types.

This function queries PIASOmarkerDB to find which cell types are associated with the input genes, ranking results by matched gene count and specificity.

Parameters

genes — list of str, pd.DataFrame, or dict

Gene input. Supports three formats:

list of str: Single list of gene symbols. Returns: pd.DataFrame with analysis results.
pd.DataFrame: Columns are clusters/cell types, rows are genes. Ideal for COSG output: pd.DataFrame(adata.uns['cosg']['names']).head(50). Returns: tuple (results_dict, top_hits_dict).
dict: {cluster_name: [gene_list]}. Returns: tuple (results_dict, top_hits_dict).

n_top_genes — int, optional

For DataFrame/dict input: only use top N genes per column/key. Default: 50. Useful for COSG results which may rank many genes.

species — str, optional

Filter results by species (e.g., “Human”, “Mouse”).

tissue — str, optional

Filter results by tissue.

studies — str or list of str, optional

Study/studies to include in analysis. Only cell types from these studies will be considered. Study names are validated against PIASOmarkerDB. Cannot overlap with exclude_studies. Default: None (use all studies).

min_genes — int, optional

Minimum number of genes that must match a cell type. Default: 1.

exclude_cell_types — list of str, optional

Cell types to exclude from results.

exclude_studies — list of str, optional

Studies to exclude from results. Cannot overlap with studies.

Returns

pd.DataFrame

For single list input: DataFrame with columns: cell_type, study_publication, species, tissue, condition, matched_gene_count, matched_genes, avg_specificity.

tuple (dict, dict) For DataFrame or dict input:

results_dict: {cluster_name: result_DataFrame}
top_hits_dict: {cluster_name: "predicted_cell_type"} (or “Unassigned” if no matches found)

Raises

ValidationError

If studies parameter contains invalid study names, or if studies and exclude_studies have overlapping values.

Examples

Single gene list (mouse cortex L6 markers):

>>> import piaso
>>> query_genes = ["Syt6", "Tle4", "Hs3st4", "Fezf2", "Foxp2", "Col12a1"]
>>> df = piaso.tl.analyzeMarkers(query_genes)
>>> print(df.head())
# Top hit: EN-L6-CT from WangKriegstein2025

With specific study filter:

>>> results, top_hits = piaso.tl.analyzeMarkers(
...     cosg_marker_df,
...     n_top_genes=50,
...     min_genes=5,
...     studies=['AllenWholeMouseBrain_isocortex'],
...     species="Mouse"
... )
>>> print(top_hits)
{'L2-3 IT': '007 L2/3 IT CTX Glut', 'PV': '052 Pvalb Gaba', ...}

Dictionary input (microglia and L6 CT markers):

>>> gene_sets = {
...     'Cluster_0': ['Cx3cr1', 'P2ry12', 'Tmem119', 'Csf1r', 'Trem2'],
...     'Cluster_1': ['Syt6', 'Tle4', 'Hs3st4', 'Fezf2', 'Foxp2'],
... }
>>> results, top_hits = piaso.tl.analyzeMarkers(gene_sets)
>>> print(top_hits)
{'Cluster_0': 'Microglia', 'Cluster_1': 'EN-L6-CT'}

COSG integration workflow:

>>> import cosg
>>> import pandas as pd
>>>
>>> # Run COSG
>>> cosg.cosg(adata, key_added='cosg', groupby='leiden')
>>>
>>> # Get top 50 markers per cluster
>>> cosg_marker_df = pd.DataFrame(adata.uns['cosg']['names']).head(50)
>>>
>>> # Analyze with PIASOmarkerDB
>>> results, top_hits = piaso.tl.analyzeMarkers(
...     cosg_marker_df,
...     n_top_genes=50,
...     species="Mouse"
... )
>>>
>>> # Add annotations to AnnData
>>> adata.obs['cell_type_predicted'] = adata.obs['leiden'].map(top_hits)

See Also

queryPIASOmarkerDB : Direct marker queries

Notes

For COSG results, the DataFrame columns are cluster/cell type names and rows contain the ranked marker genes. Only the top n_top_genes genes per column are used for analysis.

When studies is provided, study names are validated against PIASOmarkerDB. If an invalid study name is provided, a ValidationError is raised with instructions to list available studies.

PIASOmarkerDB website: https://piaso.org/piasomarkerdb/

calculateScoreParallel

calculateScoreParallel(
    adata,
    gene_set: Union[dict, list, pandas.core.frame.DataFrame],
    gene_weights=None,
    score_method: Literal['scanpy', 'piaso'] = 'piaso',
    random_seed: int = 1927,
    score_layer: Optional[str] = None,
    max_workers: Optional[int] = None,
    return_pvals: bool = False,
    precomputed_knn: numpy.ndarray = None,
    verbosity: int = 0,
    modality: str = 'RNA',
    cytome_layer: str = 'counts',
    batch_size: int = 1024,
    cell_mask=None,
    score_chunk_size: Optional[int] = None,
    max_score_chunk_bytes: int = 268435456,
    max_score_batch_cache_bytes: int = 536870912,
)

Signature defaults

adata, gene_set, gene_weights=None, score_method='piaso', random_seed=1927, score_layer=None, max_workers=None, return_pvals=False, precomputed_knn=None, verbosity=0, modality='RNA', cytome_layer='counts', batch_size=1024, cell_mask=None, score_chunk_size=None, max_score_chunk_bytes=268435456, max_score_batch_cache_bytes=536870912

Compute gene set scores in parallel using shared memory for efficiency.

This function processes multiple gene sets in parallel, computing enrichment scores for each gene set across all cells in the AnnData object. When using the ‘piaso’ scoring method, it uses a vectorized batched approach (score() with multi-set mode) that precomputes gene-level statistics once and scores all gene sets in a single pass, which is significantly faster and more memory-efficient than scoring each set independently. For the ‘scanpy’ method, it uses shared memory to pass the expression matrix to worker processes.

Parameters

adata — AnnData

The input AnnData object containing gene expression data.

gene_set — dict, list of lists, or pandas.DataFrame

A collection of gene sets to score. Supported formats:

dict: Keys are gene set names, values are lists of gene names.
list of lists: Each sublist contains gene names for one gene set. Gene sets will be named “GeneSet_0”, “GeneSet_1”, etc.
pandas.DataFrame: Each column represents a gene set, with column names as gene set names and gene names as values.

gene_weights — list of array-like, optional

One weight vector per gene set, in the same order as gene_set. Passed straight to :func:piaso.tl.score. Only valid with score_method='piaso'.

score_method — {‘scanpy’, ‘piaso’}, default ‘piaso’

The method used for gene set scoring.

‘scanpy’: Uses Scanpy’s built-in gene set scoring method.
‘piaso’: Uses the PIASO’s gene set scoring method, which is more robust to sequencing depth variations and provides p-values.

random_seed — int, default 1927

Random seed for reproducibility.

score_layer — str or None, default None

Layer of the AnnData object to use. If None, adata.X is used.

max_workers — int or None, default None

Number of parallel worker processes to use. If None, defaults to the number of CPU cores available. Only used when score_method=‘scanpy’.

return_pvals — bool, default False

Whether to return -log10(p-values) when using ‘piaso’ method. Only applicable when score_method=‘piaso’. If True, returns a third array containing p-values.

verbosity — int, default 0

Level of verbosity for progress reporting.

0: Silent (no progress bar)
0: Show progress bar during parallel computation

precomputed_knn — optional

Reuse a control-gene neighbour structure from an earlier call instead of rebuilding it.

modality — str, default ‘RNA’

Cytome only: the modality to score.

cytome_layer — str, default ‘counts’

Cytome only: the layer to read.

batch_size — int, default 1024

Cytome only: rows per streaming chunk.

cell_mask — ndarray, optional

Cytome only: boolean mask or sorted indices of the cells to score.

score_chunk_size — int, optional

Cytome only: rows handed to the scoring kernel per call; None derives it from max_score_chunk_bytes.

max_score_chunk_bytes — int, default 256 MB

Cytome only: memory budget behind score_chunk_size, total across the concurrent workers.

max_score_batch_cache_bytes — int, default 512 MB

Cytome only: budget for holding a batch’s chunks between the two streaming passes so the second does not re-read them.

Returns

score_matrix — np.ndarray

A 2D array of shape (n_cells, n_gene_sets) where each column contains the scores for one gene set across all cells.

gene_set_names — list of str

The names of the gene sets, in the same order as columns in score_matrix.

nlog10_pval_matrix — np.ndarray, optional

Only returned when score_method=‘piaso’ and return_pvals=True. A 2D array of shape (n_cells, n_gene_sets) containing -log10(p-values) for each gene set score. Returns None if p-values are not available.

Examples

>>> import anndata
>>> import numpy as np
>>> import piaso
>>>
>>> # Load example data
>>> adata = anndata.read_h5ad('pbmc3k.h5ad')
>>>
>>> # Define gene sets
>>> gene_sets = {
...     'T_cell_markers': ['CD3D', 'CD3E', 'CD8A'],
...     'B_cell_markers': ['CD79A', 'CD79B', 'MS4A1']
... }
>>>
>>> # Compute scores using Scanpy method
>>> scores, names = piaso.tl.calculateScoreParallel(
...     adata,
...     gene_set=gene_sets,
...     score_method='piaso',
...     verbosity=1
... )
>>>
>>> # Add scores to AnnData object
>>> for i, name in enumerate(names):
...     adata.obs[f'{name}_score'] = scores[:, i]

calculateScoreParallel_multiBatch

calculateScoreParallel_multiBatch(
    adata,
    batch_key: str,
    marker_gene: pandas.core.frame.DataFrame,
    marker_gene_n_groups_indices: list,
    score_method: Literal['scanpy', 'piaso'],
    score_layer: str = None,
    max_workers: int = 8,
    n_concurrent_batches: int = None,
    random_seed: int = 1927,
)

Signature defaults

adata, batch_key, marker_gene, marker_gene_n_groups_indices, score_method, score_layer=None, max_workers=8, n_concurrent_batches=None, random_seed=1927

Calculate gene set scores for each adata batch in parallel using shared memory. Different marker gene sets will be calculated in parallel as well.

Parameters

adata — AnnData

Annotated data matrix.

batch_key — str

The key in adata.obs used to identify batches.

marker_gene — DataFrame

The marker gene DataFrame.

marker_gene_n_groups_indices — list

Indices specifying the marker gene set group boundaries, used for score normalization within each marker gene set group.

max_workers — int

Maximum number of parallel workers to use (total threads).

score_layer — str

The layer of adata to use for scoring.

score_method — {‘scanpy’, ‘piaso’}, optional

The method used for gene set scoring. Must be either ‘scanpy’ (default) or ‘piaso’.

‘scanpy’: Uses the Scanpy’s built-in gene set scoring method.
‘piaso’: Uses the PIASO’s gene set scoring method, which is more robust to sequencing depth variations.

n_concurrent_batches — int, optional

Number of batches to process concurrently via ThreadPoolExecutor. If None, auto-determined based on max_workers and number of batches. Only used when score_method=‘piaso’. Default is None.

random_seed — int, optional

Random seed for reproducibility. Default is 1927.

Returns

tuple

list: A list of normalized score arrays for each batch.
list: A list of cell barcodes for each batch.
list: A list of gene set names.

Examples

>>> import anndata
>>> import piaso
>>> adata = anndata.read_h5ad('example_data.h5ad')
>>> score_list, cellbarcode_info, gene_set_names = piaso.tl.calculateScoreParallel_multiBatch(
...     adata=adata,
...     batch_key='batch',
...     marker_gene=marker_gene,
...     marker_gene_n_groups_indices=marker_gene_n_groups_indices,
...     score_layer='piaso',
...     max_workers=8
... )
>>> print(score_list)
>>> print(cellbarcode_info)

compute_tfidf_stats

compute_tfidf_stats(
    source,
    measurement: str = 'counts',
    batch_size: int = 1024,
    scale_factor: float = 10000.0,
    modality: str = 'ATAC',
    write_to_metadata: bool = True,
)

Signature defaults

source, measurement='counts', batch_size=1024, scale_factor=10000.0, modality='ATAC', write_to_metadata=True

Compute TF-IDF statistics without materializing the full matrix.

Pass 1 only: accumulates cell_depth and peak_depth, computes idf. The returned arrays can be passed to runSVD via tfidf_params for inline TF-IDF application during SVD — eliminating the persistent TF-IDF layer entirely.

Parameters

source — CytomeDataset or str

Cytome dataset or path.

measurement — str

Input measurement name (default ‘counts’).

batch_size — int

Cells per chunk.

scale_factor — float

TF-IDF scale factor.

modality — str

Modality prefix (default ‘ATAC’).

write_to_metadata — bool, default True

Cache the statistics in the cytome’s metadata under '{modality}_tfidf_params' so later calls reuse them.

Returns

dict

{‘cell_depth’: ndarray, ‘idf’: ndarray, ‘scale_factor’: float}

cospecificity_trans

_cospecificity_trans_impl(
    data,
    groupby: 'str' = 'leiden',
    pairs: 'dict | None' = None,
    metric: 'str' = 'weighted_cosine',
    modality: 'str' = 'RNA',
    cosg_mu: 'float' = 1.0,
    cosg_expressed_pct: 'float' = 0.1,
    cosg_layer: 'str' = 'counts',
    cosg_batch_size=None,
    use_cached: 'bool' = True,
    cosg_cache_key: 'str | None' = None,
    specificity: 'pd.DataFrame | None' = None,
    min_sign: 'float' = 0.0,
    verbose: 'int' = 1,
)

Signature defaults

data, groupby='leiden', pairs=None, metric='weighted_cosine', modality='RNA', cosg_mu=1.0, cosg_expressed_pct=0.1, cosg_layer='counts', cosg_batch_size=None, use_cached=True, cosg_cache_key=None, specificity=None, min_sign=0.0, verbose=1

Trans co-specificity between candidate feature pairs (e.g. TF→target).

Parameters

data

AnnData or cytome.Dataset (RNA).

pairs

{source: [target, ...]} candidate edges (e.g. the cistrome’s motif-supported TF→gene map). Only these edges are scored (keeps it tractable). If None, all source×target pairs over the features are used (expensive — intended for small inputs/tests).

metric

Per-cell-type pairwise metric (weighted_cosine default; also geomean / cosine / outer) — reuses :func:_pairwise_metric.

specificity

Optional precomputed S [features × cell_types]; skips COSG.

min_sign

Keep only edges whose specificity-profile Pearson correlation (activation sign) > min_sign. Default 0 → positive (activating) edges only (design Q6/Q7: no repression inference).

Returns

dict

{"edges": DataFrame[source,target,cosine,sign,cospec_max,best_celltype], "celltypes": [...], "per_celltype": ndarray[n_edges, n_celltypes], "S": DataFrame}. per_celltype[e, t] is the metric for edge e in cell type t — the basis for cell-type-specific regulons.

ensure_infog_params

_ensure_infog_params(ds, modality, use_cached_stats=True, batch_size=2048)

Signature defaults

ds, modality, use_cached_stats=True, batch_size=2048

Return cached {modality}_infog_params (or compute now and cache). Falls through to legacy infog_params once with a DeprecationWarning when only the unprefixed key is present. A cached payload whose cell_depth length no longer matches ds.n_cells is treated as a miss (self-heals a cytome filtered before the cached-stats invalidation fix).

ensure_tfidf_params

_ensure_tfidf_params(ds, modality, use_cached_stats=True, batch_size=2048)

Signature defaults

ds, modality, use_cached_stats=True, batch_size=2048

Return cached {modality}_tfidf_params (or compute now and cache).

Delegates to _runTFIDF._load_or_compute_tfidf_stats so the runSVD auto_tfidf=True path, COSG layer='tfidf', and plotting cytome_layer='tfidf' all share the same cache-or-compute helper.

getMarkers

queryPIASOmarkerDB(
    gene: 'Union[str, List[str]]' = None,
    cell_type: 'Union[str, List[str]]' = None,
    study: 'str' = None,
    species: 'str' = None,
    tissue: 'str' = None,
    condition: 'str' = None,
    min_score: 'float' = None,
    max_score: 'float' = None,
    limit: 'int' = None,
    as_dict: 'Union[bool, str]' = False,
    list_studies: 'bool' = False,
    list_cell_types: 'bool' = False,
    list_genes: 'bool' = False,
)

Signature defaults

gene=None, cell_type=None, study=None, species=None, tissue=None, condition=None, min_score=None, max_score=None, limit=None, as_dict=False, list_studies=False, list_cell_types=False, list_genes=False

Query PIASOmarkerDB for cell type marker genes.

This is the main entry point for accessing PIASOmarkerDB through PIASO.

Parameters

gene — str or list of str, optional

Gene symbol(s) to filter by.

cell_type — str or list of str, optional

Cell type(s) to filter by.

study — str, optional

Study/publication to filter by.

species — str, optional

Species to filter by (e.g., “Human”, “Mouse”).

tissue — str, optional

Tissue to filter by.

condition — str, optional

Condition to filter by.

min_score — float, optional

Minimum specificity score (>= 0).

max_score — float, optional

Maximum specificity score (>= 0).

limit — int, optional

Maximum results to return. Default: None (no limit).

as_dict — bool or {‘both’}, optional

True returns the {cell_type: [genes]} dictionary alone — the form predictCellTypeByMarker takes. False (the default) returns the DataFrame. 'both' returns the old (DataFrame, dict) pair and is deprecated; it will be removed in the next release.

list_studies — bool, optional

If True, return list of available study names instead of markers. Default: False.

list_cell_types — bool, optional

If True, return list of available cell types instead of markers. Default: False.

list_genes — bool, optional

If True, return list of unique gene symbols instead of markers. Default: False.

Returns

pd.DataFrame

Marker query results (default).

dict

If as_dict=True: {cell_type: [genes]}.

tuple (pd.DataFrame, dict) If as_dict=‘both’ (deprecated).

list of str If list_studies=True, list_cell_types=True, or list_genes=True.

Examples

Query marker genes:

>>> import piaso
>>> df = piaso.tl.queryPIASOmarkerDB(gene="Foxp2", species="Mouse")
>>> df = piaso.tl.queryPIASOmarkerDB(gene=["Foxp2", "Syt6", "Tle4"])

Get the marker dictionary, ready for predictCellTypeByMarker:

>>> marker_dict = piaso.tl.queryPIASOmarkerDB(
...     study="AllenWholeMouseBrain_isocortex",
...     species="Mouse",
...     as_dict=True
... )
>>> print(f"Cell types in dict: {len(marker_dict)}")
>>> piaso.tl.predictCellTypeByMarker(adata, marker_gene_set=marker_dict)

List available studies:

>>> studies = piaso.tl.queryPIASOmarkerDB(list_studies=True)
>>> print(f"Total studies: {len(studies)}")

List cell types:

>>> cell_types = piaso.tl.queryPIASOmarkerDB(list_cell_types=True, species="Mouse")

See Also

analyzeMarkers : Analyze gene lists for cell type inference PIASOmarkerDB : Client class for advanced usage

Notes

PIASOmarkerDB website: https://piaso.org/piasomarkerdb/

infog

infog(
    data=<unset>,
    copy: bool = False,
    inplace: bool = False,
    n_top_genes: int = 3000,
    key_added: str = 'infog',
    key_added_highly_variable_gene: str = 'highly_variable',
    trim: bool = True,
    verbosity: int = 1,
    layer: Optional[str] = None,
    streaming: bool = False,
    batch_size: int = 1024,
    save_layer: bool = False,
    dtype: Optional[str] = None,
    modality: str = 'RNA',
    return_info: bool = False,
    allow_non_integer: bool = False,
    source=<unset>,
    adata=<unset>,
)

Signature defaults

data=<unset>, copy=False, inplace=False, n_top_genes=3000, key_added='infog', key_added_highly_variable_gene='highly_variable', trim=True, verbosity=1, layer=None, streaming=False, batch_size=1024, save_layer=False, dtype=None, modality='RNA', return_info=False, allow_non_integer=False, source=<unset>, adata=<unset>

INFOG normalization of single-cell RNA sequencing data.

Supports three modes:

infog(adata) — standard in-memory (existing behavior, unchanged)
infog(adata, streaming=True) — streaming from in-memory AnnData chunks
infog(“path.cytome”) — streaming from on-disk cytome dataset
infog(cytome_dataset) — streaming from already-opened cytome object

Parameters

source — AnnData, CytomeDataset, or str

AnnData object, cytome Dataset object, or path to .cytome file.

modality — str, default ‘RNA’

Cytome modality to read. Only meaningful when source is a cytome path or Dataset; ignored for AnnData inputs. Use ‘GA’ to compute INFOG on gene-activity matrices (e.g. after piaso.tl.inferGeneActivity). The math is modality-agnostic; passing ‘ATAC’ / ‘tiles’ is technically supported but biologically unusual.

streaming — bool, default=False

If True and source is AnnData, use streaming mode. Ignored if source is str or CytomeDataset (these always use streaming).

batch_size — int, default=1024

Number of cells per chunk in streaming mode. Ignored in standard mode.

save_layer — bool, default=False

If True and source is cytome, write the full INFOG-normalized matrix to the cytome file (layer key_added). Default (False) is lazy mode: only normalization parameters are saved, and normalization is applied on-the-fly during downstream operations like SVD.

allow_non_integer — bool, default False

By default INFOG refuses input whose values are not integers, because its dispersion model is defined on raw UMI counts and silently returns meaningless numbers on normalized, log-transformed or scaled data. Set True to run anyway — appropriate for Smart-seq2 TPM/FPKM, imputed or already-corrected matrices. Ignored when layer/infog_layer is given, since naming a layer already answers the question. [all other parameters unchanged from original infog()]

data — AnnData, CytomeDataset or str

The data: an AnnData, an open cytome, or a path to one. The preferred name for source.

copy — bool, default False

AnnData: return a normalised copy instead of editing in place.

inplace — bool, default False

AnnData: store the normalised matrix in adata.X rather than in adata.layers[key_added].

n_top_genes — int, default 3000

Number of highly variable genes to select.

key_added — str, default ‘infog’

Name of the normalised layer.

key_added_highly_variable_gene — str, default ‘highly_variable’

Name of the boolean gene-selection column.

trim — bool, default True

Trim extreme normalised values, as INFOG does by default.

verbosity — int, default 1

0 is silent.

layer — str, optional

AnnData: the layer holding raw counts; None is adata.X.

dtype — str, optional

Storage dtype of the normalised layer; None keeps the package default (float32).

return_info — bool, default False

Also return a dict with the normalisation parameters and the selected gene indices, so a caller can normalise more data the same way.

Other Parameters

source, adata Deprecated aliases for data; passing either warns.

infog_svd

infog_svd(
    source,
    copy: bool = False,
    n_components: int = 50,
    use_highly_variable: bool = True,
    n_top_genes: int = 3000,
    verbosity: int = 0,
    batch_key: Optional[str] = None,
    random_state: Optional[int] = 1927,
    scale_data: bool = False,
    n_iter: int = 7,
    infog_trim: bool = True,
    key_added: str = 'X_svd',
    layer: Optional[str] = 'infog',
    infog_layer: Optional[str] = None,
    streaming: bool = False,
    batch_size: int = 1024,
    allow_non_integer: bool = False,
)

Signature defaults

source, copy=False, n_components=50, use_highly_variable=True, n_top_genes=3000, verbosity=0, batch_key=None, random_state=1927, scale_data=False, n_iter=7, infog_trim=True, key_added='X_svd', layer='infog', infog_layer=None, streaming=False, batch_size=1024, allow_non_integer=False

INFOG normalization → HVG selection → SVD in one call.

Performs INFOG normalization, selects highly variable genes, and runs truncated SVD for dimensionality reduction in a single function call.

Supports three modes:

infog_svd(adata) — standard in-memory (existing behavior)
infog_svd(adata, streaming=True) — streaming from in-memory AnnData
infog_svd(“path.cytome”) — streaming from on-disk cytome dataset
infog_svd(cytome_dataset) — streaming from already-opened cytome object allow_non_integer : bool, default False By default INFOG refuses input whose values are not integers, because its dispersion model is defined on raw UMI counts and silently returns meaningless numbers on normalized, log-transformed or scaled data. Set True to run anyway — appropriate for Smart-seq2 TPM/FPKM, imputed or already-corrected matrices. Ignored when layer/infog_layer is given, since naming a layer already answers the question.

Parameters

source — AnnData, CytomeDataset or str

The data: an AnnData, an open cytome, or a path to one.

copy — bool, default False

AnnData: return a copy with the embedding instead of editing in place.

n_components — int, default 50

SVD dimensions.

use_highly_variable — bool, default True

Restrict the SVD to the selected genes.

n_top_genes — int, default 3000

Genes INFOG selects.

verbosity — int, default 0

1 reports each step.

batch_key — str, optional

Kept for API compatibility; the selection no longer depends on it.

random_state — int, default 1927

Seed for the randomised SVD.

scale_data — bool, default False

Scale features before the SVD (AnnData path only).

n_iter — int, default 7

Rounds of the randomised solver.

infog_trim — bool, default True

Passed to INFOG as trim.

key_added — str, default ‘X_svd’

Name of the embedding.

layer — str, default ‘infog’

'infog' (the default) normalises with INFOG first, reading the raw counts from infog_layer (None: adata.X), and selects genes as INFOG does. Any other value names a layer to factorise as it is, with the object’s highly_variable column if it has one and every gene otherwise; None means adata.X as it is.

infog_layer — str, optional

The layer holding raw counts for INFOG (None: adata.X). Used only when layer='infog'.

streaming — bool, default False

AnnData: use the chunked streaming path; a cytome always streams.

batch_size — int, default 1024

Rows per chunk on the streaming path.

allow_non_integer — bool, default False

Let INFOG accept a matrix that does not hold integer counts.

leiden

leiden(
    data,
    resolution=1.0,
    modality='RNA',
    n_iterations=10,
    random_state=42,
    key_added='leiden',
    adjacency_key=None,
    neighbors_key=None,
    knn_result=None,
    cell_mask=None,
    backend='auto',
    n_threads=0,
)

Signature defaults

data, resolution=1.0, modality='RNA', n_iterations=10, random_state=42, key_added='leiden', adjacency_key=None, neighbors_key=None, knn_result=None, cell_mask=None, backend='auto', n_threads=0

Leiden clustering of the neighbour graph.

Reads the connectivities matrix and runs the Leiden algorithm: PIASO’s own implementation, in its Rust extension, or igraph’s community_leiden().

Parameters

data — AnnData, cytome.Dataset, or the path of a cytome

If AnnData: reads from obsp, stores in obs. If cytome.Dataset or a path: reads connectivities from cytome graphs, stores in cells.

resolution — float

Resolution parameter controlling cluster granularity: higher gives more and smaller clusters.

n_iterations — int, default 10

How many times the whole algorithm runs, each run starting from the partition the one before it ended with. A negative number repeats until a run changes nothing, which with random refinement can take long. More runs raise the modularity a little and make two seeds agree a little more; on 2.3 million cells 2 runs give 96 to 99 % of the partition of 10 in under half the time.

random_state — int or None, default 42

The seed. With backend='rust' one graph and one seed give one partition, on any machine and any number of threads; None refines without a random choice. With backend='igraph' it seeds igraph’s generator, and repeated calls on one machine agree.

key_added — str

Column name to store cluster labels.

neighbors_key — str, optional

Prefix of the neighbors graph to use (matches the key_added passed to piaso.tl.neighbors). neighbors_key='SVD' reads 'SVD_connectivities'; None / 'neighbors' reads the un-prefixed 'connectivities'. Mirrors umap’s neighbors_key.

adjacency_key — str, optional

Full graph name to read (escape hatch / overrides neighbors_key). Defaults to the connectivities graph resolved from neighbors_key.

knn_result — dict, optional

Result dict from neighbors() with ‘connectivities’. Used for the in-memory ndarray / cell_mask path to avoid re-reading from disk.

modality — str, default ‘RNA’

Cytome only: the modality the graph and the labels belong to.

cell_mask — ndarray, optional

Cytome only: boolean mask or sorted indices of the cells the graph was built on with the same mask in neighbors.

backend — {‘auto’, ‘rust’, ‘igraph’}, default ‘auto’

'auto' is 'rust' when the extension provides it, else 'igraph'. 'rust': PIASO’s Leiden, parallel and deterministic. One graph and one random_state give one partition: the number of threads, the cores of the machine and the layout of the graph in memory do not enter the result. Same objective and the same randomised refinement as igraph, and the same modularity; random_state=None refines greedily, with no random choice. It releases Python’s lock while it runs. The graph must be symmetric, and a directed one is refused. The graph passed in is not modified. Ctrl-C stops a run.

Everything in PIASO that clusters calls this function, so leiden_local and runGDR follow the same backend.

Changed in 1.2.6: 'rust' is the default, and the partitions differ from those of 1.2.5 by about as much as another seed would change them. Pass backend='igraph' for the partitions of earlier versions.

n_threads — int, default 0

'rust' only: threads to use. 0 uses all the cores this process may use (on a cluster: the job’s allocation), or the number in the environment variable RAYON_NUM_THREADS if it is set. The labels are the same for every value. A graph of fewer than 32,768 nodes is clustered on one thread whatever is asked, which is faster there.

Returns

None, or np.ndarray when the labels cannot be written Nothing, on both backends: the labels are written to adata.obs[key_added] or to ds.cells[key_added] and read back from there. The two paths that cannot write still return the label array — cell_mask (the labels cover the masked cells only, so writing them full-length would mis-align) and data=None (no object to write to).

Changed in 1.2.4: the AnnData path used to return the labels as well as writing them, so the same call had two different contracts depending on the backend and echoed an array of strings in a notebook. Read them from adata.obs[key_added].

leiden_local

leiden_local(
    adata,
    clustering_type: str = 'each',
    groupby: str = 'Leiden',
    groups: Optional[Sequence[str]] = None,
    resolution: float = 0.25,
    batch_key: Optional[Sequence[str]] = None,
    key_added: str = 'Leiden_local',
    dr_method: str = 'X_svd',
    gdr_resolution: float = 1.0,
    copy: bool = False,
    modality: Optional[str] = None,
    cytome_layer: str = 'counts',
    max_nnz_percentile: float = 20.0,
    n_components: int = 30,
    svd_method: str = 'auto',
    n_iter: Optional[int] = None,
    svd_tol: Optional[float] = None,
    n_neighbors: int = 15,
    batch_size: int = 1024,
    random_state: int = 10,
    min_cells: int = 20,
    verbosity: int = 1,
    n_threads: int = 0,
    backend: str = 'auto',
)

Signature defaults

adata, clustering_type='each', groupby='Leiden', groups=None, resolution=0.25, batch_key=None, key_added='Leiden_local', dr_method='X_svd', gdr_resolution=1.0, copy=False, modality=None, cytome_layer='counts', max_nnz_percentile=20.0, n_components=30, svd_method='auto', n_iter=None, svd_tol=None, n_neighbors=15, batch_size=1024, random_state=10, min_cells=20, verbosity=1, n_threads=0, backend='auto'

Leiden clustering within each group of a coarse clustering.

Every group is re-embedded in its own low-dimensional space and clustered there, so structure the whole-dataset embedding merged is recovered. Works on an AnnData in memory and on a cytome (open Dataset or path) streaming, per group, in bounded memory.

Parameters

adata — AnnData, cytome.Dataset or str

The data; a cytome may be given as a path.

clustering_type — {‘each’, ‘all’}, default ‘each’

'each' clusters every selected group on its own; 'all' pools the selected groups and clusters them together (sub-labels prefixed M-). A cytome supports 'each' only.

groupby — str, default ‘Leiden’

The coarse labels: an obs column, or a cells column on a cytome.

groups — sequence of str, optional

Which groups to cluster; all of them when None.

resolution — float, default 0.25

Leiden resolution inside each group; higher gives more clusters.

batch_key — str, optional

Batch column for the Harmony variants of dr_method.

key_added — str, default ‘Leiden_local’

Where the labels go: obs[key_added] or cells[key_added], as <group>-<k>.

dr_method — str, default ‘X_svd’

How each group is embedded:

'X_svd': runSVD on the group’s selected features. On an ATAC or tiles cytome the selection is made per group by cumulative non-zeros (max_nnz_percentile) with TF-IDF; on an RNA or GA cytome and on an AnnData it is the highly_variable column, or every feature if there is none (said once). An ATAC AnnData should be TF-IDF-normalised first (run_TFIDF).
'X_svd_full': runSVD on every feature with at least one count, values as stored.
'X_pca': INFOG normalisation and selection, then SVD.
'X_gdr', 'X_pca_harmony', 'X_svd_full_harmony': as their names say; Harmony needs batch_key.

On a cytome, 'X_svd' and 'X_svd_full' stream per group; the other methods convert the store to an AnnData first.

gdr_resolution — float, default 1.0

Leiden resolution of the temporary clustering GDR is built on (dr_method='X_gdr').

copy — bool, default False

AnnData only: return a modified copy instead of editing in place.

modality — str, optional

Cytome only: which modality to cluster ('ATAC', 'tiles', 'RNA', 'GA'). Defaults to 'ATAC' on the streaming path and 'RNA' for the AnnData conversion.

cytome_layer — str, default ‘counts’

Cytome only: the layer to read. For RNA that is usually a materialised 'infog'.

max_nnz_percentile — float, default 20.0

ATAC/tiles streaming only: the cumulative-non-zero percentile that selects a group’s peaks.

n_components — int, default 30

SVD dimensions per group, capped at one less than the group’s cell count and than its selected feature count (and, for the AnnData path, at the square root of the cell count, as before).

svd_method — {‘auto’, ‘power’, ‘krylov’}, default ‘auto’

The per-group solver; see runSVD.

n_iter — int, optional

Rounds for 'power'/'krylov', a ceiling for 'auto'; see runSVD.

svd_tol — float, optional

Convergence tolerance passed to runSVD.

n_neighbors — int, default 15

Neighbours for the per-group graph (capped below the group size).

batch_size — int, default 1024

Rows per chunk for the Python reader, where it is used.

random_state — int, default 10

Seed for the SVD, the graph and Leiden.

min_cells — int, default 20

Cytome streaming only: a group smaller than this is kept as one sub-cluster <group>-0 rather than re-embedded.

verbosity — int, default 1

1 reports each group’s outcome.

n_threads — int, default 0

Threads for each group’s SVD and Leiden; 0 uses all the cores this process may use. The labels do not depend on it.

backend — {‘auto’, ‘rust’, ‘igraph’}, default ‘auto’

The Leiden of every group, as in :func:piaso.tl.leiden.

Returns

AnnData or None The copy when copy=True; otherwise None after writing the labels in place.

Notes

On the cytome streaming path a group whose chain fails is kept as one sub-cluster <group>-0 with a RuntimeWarning naming the error, so one bad group does not lose the run; a result full of -0 labels is the sign to read the warnings.

neighbors

neighbors(
    data,
    use_rep='X_svd',
    n_neighbors=15,
    modality='RNA',
    metric='euclidean',
    random_state=42,
    key_added=None,
    cell_mask=None,
)

Signature defaults

data, use_rep='X_svd', n_neighbors=15, modality='RNA', metric='euclidean', random_state=42, key_added=None, cell_mask=None

Build kNN graph and compute fuzzy simplicial set connectivities.

Uses pynndescent for approximate nearest neighbor search and umap’s fuzzy_simplicial_set for UMAP-compatible connectivities.

Parameters

data — AnnData or cytome.Dataset

If AnnData: reads from obsm, stores in obsp/uns. If cytome.Dataset: reads from embeddings, stores graphs in cytome.

use_rep — str

Embedding name. For AnnData: key in obsm. For cytome: embedding name.

n_neighbors — int

Number of nearest neighbors.

metric — str

Distance metric for pynndescent.

random_state — int

Random seed for reproducibility.

key_added — str, optional

Prefix for the stored graph names. None (default) writes the un-prefixed connectivities / distances (+ n_neighbors metadata); 'SVD' writes SVD_connectivities / SVD_distances. The legacy value 'neighbors' is a back-compat alias for “no prefix”. Pass the same string as neighbors_key to leiden / umap.

modality — str, default ‘RNA’

Cytome only: the modality use_rep was stored under, so that 'X_svd' resolves to ATAC_svd for an ATAC run.

cell_mask — ndarray, optional

Cytome only: boolean mask or sorted indices of the cells to build the graph on; the result covers those cells only.

Returns

dict or None For AnnData and the in-memory ndarray / cell_mask paths, returns a dict with ‘knn_indices’, ‘knn_dists’, ‘connectivities’, ‘distances’. For a cytome.Dataset (no cell_mask) returns None — the graph is persisted on the cytome (connectivities / distances graphs + an n_neighbors metadata entry), and piaso.tl.umap / piaso.tl.leiden read it back from there. The function is self-contained: no value passing required.

PIASOmarkerDB

PIASOmarkerDB(
    base_url: 'str' = None,
    timeout: 'int' = None,
    cache_dir: 'str | Path | None' = None,
)

Signature defaults

base_url=None, timeout=None, cache_dir=None

Python client for accessing PIASOmarkerDB.

PIASOmarkerDB is a comprehensive database of cell type marker genes with specificity scores across various tissues, species, studies, and conditions, powered by PIASO (Precise Integrative Analysis of Single-cell Omics) methodologies.

Parameters

base_url — str, optional

Base URL for the PIASOmarkerDB API. Default: “https://piaso.org/piasomarkerdb”

timeout — int, optional

Request timeout in seconds. Default: 30

cache_dir — str or Path, optional

Directory for caching downloaded markers. Default: ~/.piaso/markers

Examples

Basic usage:

>>> from piaso.tools import PIASOmarkerDB
>>> client = PIASOmarkerDB()
>>>
>>> # Query markers
>>> df = client.getMarkers(gene="Foxp2")
>>>
>>> # Get as dict
>>> marker_dict = client.getMarkers(
...     study="AllenWholeMouseBrain_isocortex",
...     as_dict=True
... )
>>>
>>> # Download markers
>>> client.downloadMarkers("markers.csv", species="Mouse")

See Also

piaso.tl.queryPIASOmarkerDB : Main query function piaso.tl.analyzeMarkers : Gene list analysis

Notes

PIASOmarkerDB provides marker genes with specificity scores computed using PIASO’s standardized methodology across multiple single-cell RNA-seq studies.

Website: https://piaso.org/piasomarkerdb/

predictCellTypeByGDR

predictCellTypeByGDR(
    adata,
    adata_ref,
    layer: str = 'log1p',
    layer_reference: str = 'log1p',
    reference_groupby: str = 'CellTypes',
    query_groupby: str = 'Leiden',
    mu: float = 10.0,
    n_genes: int = 15,
    return_integration: bool = False,
    use_highly_variable: bool = True,
    n_highly_variable_genes: int = 5000,
    n_svd_dims: int = 50,
    resolution: float = 1.0,
    scoring_method: str = None,
    key_added: str = None,
    verbosity: int = 0,
    modality: Optional[str] = None,
    cytome_layer: str = 'counts',
)

Signature defaults

adata, adata_ref, layer='log1p', layer_reference='log1p', reference_groupby='CellTypes', query_groupby='Leiden', mu=10.0, n_genes=15, return_integration=False, use_highly_variable=True, n_highly_variable_genes=5000, n_svd_dims=50, resolution=1.0, scoring_method=None, key_added=None, verbosity=0, modality=None, cytome_layer='counts'

Predicts cell types in a query dataset (adata) using the GDR dimensionality reduction method based on a reference dataset (adata_ref). To use GDR for dimensionality reduction, please refer to piaso.tl.runGDR or piaso.tl.runGDRParallel.

Parameters

adata — AnnData

The query single-cell AnnData object for which cell types are to be predicted.

adata_ref — AnnData

The reference single-cell AnnData object with known cell type annotations.

layer — str, optional (default: ‘log1p’)

The layer in adata to use for gene expression data. If None, uses the .X matrix.

layer_reference — str, optional (default: ‘log1p’)

The layer in adata_ref to use for reference gene expression data. If None, uses the .X matrix.

reference_groupby — str, optional (default: ‘CellTypes’)

The column in adata_ref.obs used to define reference cell type groupings.

query_groupby — str, optional (default: ‘Leiden’)

The column in adata.obs used to for GDR dimensionality reduction, such as clusters identified using Leiden or Louvain algorithms.

mu — float, optional (default: 10.0)

A regularization parameter for controlling the gene expression specificity, used in COSG (marker gene identification) and GDR.

n_genes — int, optional (default: 15)

The number of top specific genes per group, used in COSG and GDR.

return_integration — bool, optional (default: False)

If True, the function will return the integrated low-dimensional cell embeddings of the query dataset and reference dataset.

use_highly_variable — bool, optional (default: True)

Whether to use highly variable genes, used in GDR.

n_highly_variable_genes — int, optional (default: 5000)

The number of highly variable genes to select, if use_highly_variable is True, used in GDR.

n_svd_dims — int, optional (default: 50)

The number of dimensions to retain during SVD, used in GDR.

resolution — float, optional (default: 1.0)

Resolution parameter for clustering, used in GDR.

scoring_method — str, optional (default: None)

The method used for gene set scoring, used in GDR.

key_added — str, optional (default: None)

A key to add the predicted cell types or integration results to adata.obs. If None, CellTypes_gdr will be used.

verbosity — int, optional (default: 0)

The level of logging output. Higher values produce more detailed logs for debugging and monitoring progress.

modality — str, optional

Cytome only: the modality of the query to read; None is 'RNA'.

cytome_layer — str, default ‘counts’

Cytome only: the layer to read from the query.

Returns

None or AnnData If return_integration is True, returns an AnnData object of merged reference and query datasets with integrated cell embeddings and predicted cell types. Otherwise, updates adata in place with the predicted cell types.

Example

>>> import anndata
>>> # Load query dataset
>>> adata = anndata.read_h5ad("query_data.h5ad")
>>>
>>> # Load reference dataset with known cell type annotations
>>> adata_ref = anndata.read_h5ad("reference_data.h5ad")
>>>
>>> # Predict cell types for the query dataset
>>> piaso.tl.predictCellTypeByGDR(
>>>     adata=adata,
>>>     adata_ref=adata_ref,
>>>     layer='log1p',
>>>     layer_reference='log1p',
>>>     reference_groupby='CellTypes',
>>>     query_groupby='Leiden',
>>>     mu=10.0,
>>>     n_genes=20,
>>>     return_integration=False,
>>>     use_highly_variable=True,
>>>     n_highly_variable_genes=3000,
>>>     n_svd_dims=50,
>>>     resolution=0.8,
>>>     key_added='CellTypes_gdr',
>>>     verbosity=0
>>> )
>>>
>>> # Access the predicted cell types in the query dataset
>>> print(adata.obs['CellTypes_gdr'])

predictCellTypeByMarker

predictCellTypeByMarker(
    adata,
    marker_gene_set: Union[List, Dict, pandas.core.frame.DataFrame],
    marker_gene_weights: Optional[List, Dict, pandas.core.frame.DataFrame] = None,
    score_method: Literal['scanpy', 'piaso'] = 'piaso',
    score_layer: Optional[str] = 'infog',
    use_score: bool = True,
    max_workers: Optional[int] = None,
    smooth_prediction: bool = True,
    use_rep: str = 'X_gdr',
    k_nearest_neighbors: int = 7,
    return_confidence: bool = True,
    use_existing_adjacency_graph: bool = False,
    use_faiss: bool = False,
    key_added: str = 'CellTypes_predicted',
    extract_cell_type: bool = False,
    delimiter_cell_type: str = '-',
    inplace: bool = True,
    random_seed: int = 1927,
    verbosity: int = 1,
    n_jobs: int = -1,
    modality: Optional[str] = None,
    cytome_layer: str = 'counts',
)

Signature defaults

adata, marker_gene_set, marker_gene_weights=None, score_method='piaso', score_layer='infog', use_score=True, max_workers=None, smooth_prediction=True, use_rep='X_gdr', k_nearest_neighbors=7, return_confidence=True, use_existing_adjacency_graph=False, use_faiss=False, key_added='CellTypes_predicted', extract_cell_type=False, delimiter_cell_type='-', inplace=True, random_seed=1927, verbosity=1, n_jobs=-1, modality=None, cytome_layer='counts'

Predict cell types using marker genes and optionally smooth predictions.

This function performs cell type prediction using marker genes in two steps:

Calculate gene set scores for marker genes
Optionally smooth the predictions using k-nearest neighbors

Parameters

adata — AnnData

AnnData object containing single-cell data

marker_gene_set — list, dict, or pandas.DataFrame

Collection of marker genes for different cell types (pre-filtered/prepared)

marker_gene_weights — list, dict, or pandas.DataFrame, optional

Per-gene weights, one vector per cell type, aligned with marker_gene_set gene for gene. A dict keyed by cell type is the safest form; a DataFrame is read column by column; a plain sequence is taken in the order the sets appear. Requires score_method='piaso'.

Use this to let a strong, specific marker count for more than a weak one — COSG scores make a natural weight. The weights are used as given: there is no per-cell-type normalisation, because the score is divided by median(weights) * n_genes and its control sets carry the same weights, so rescaling one cell type’s weights cannot change which cell type wins. What matters is the weights’ shape within a set, not their size. A set whose median weight is zero is rejected rather than scored.

score_method — {‘scanpy’, ‘piaso’}, default=‘piaso’

Method to use for scoring marker genes

score_layer — str or None, default=‘infog’

Layer of the AnnData object to use for scoring

use_score — bool, default=True

Whether to use scores (True) or p-values (False) for cell type prediction

max_workers — int or None, default=None

Number of parallel workers for score calculation

smooth_prediction — bool, default=True

Whether to smooth predictions using k-nearest neighbors

use_rep — str, default=‘X_gdr’

Key in adata.obsm containing the low-dimensional embedding to use for neighbor search

k_nearest_neighbors — int, default=7

Number of neighbors to consider for smoothing

return_confidence — bool, default=True

Whether to return confidence scores for smoothed predictions

use_existing_adjacency_graph — bool, default=False

Whether to use existing neighborhood graph if available

use_faiss — bool, default=False

Whether to use FAISS for faster neighbor search

key_added — str, default=‘CellTypes_predicted’

Key to use for storing cell type predictions in adata.obs

extract_cell_type — bool, default=False

Whether to extract cell type name by removing suffix after delimiter

delimiter_cell_type — str, default=’-’

Delimiter to use when extracting cell type names (only used if extract_cell_type=True)

inplace — bool, default=True

Whether to modify adata in place or return a copy

random_seed — int, default=1927

Random seed for reproducibility

verbosity — int, default=1

Level of verbosity (0=quiet, 1=basic info, 2=detailed)

n_jobs — int, default=-1

Number of jobs for parallel processing during smoothing

modality — str, optional

Cytome only: the modality of the query to read; None is 'RNA'.

cytome_layer — str, default ‘counts’

Cytome only: the layer to read from the query.

Returns

If inplace=False: AnnData: Copy of adata with cell type predictions added If inplace=True: None, but adata is modified in place

Examples

>>> import piaso
>>>
>>> # Basic usage
>>> piaso.tl.predictCellTypeByMarker(
...     adata,
...     marker_gene_set=cosgMarkerDB,
...     score_method='piaso',
...     use_score=False,
...     smooth_prediction=True,
...     inplace=True
... )

projectGDR

projectGDR(
    data,
    reference,
    mode: Literal['reference', 'self'] = 'reference',
    layer: Optional[str] = None,
    counts_layer: Optional[str] = None,
    key_added: str = 'X_gdr',
    modality: str = 'RNA',
    reference_modality: Optional[str] = None,
    max_workers: int = 8,
    batch_size: int = 1024,
    min_cells_self_mode: int = 500,
    novelty_k: int = 15,
    novelty_quantile: float = 0.99,
    write_to_cytome: bool = True,
    copy: bool = False,
    verbosity: int = 1,
)

Signature defaults

data, reference, mode='reference', layer=None, counts_layer=None, key_added='X_gdr', modality='RNA', reference_modality=None, max_workers=8, batch_size=1024, min_cells_self_mode=500, novelty_k=15, novelty_quantile=0.99, write_to_cytome=True, copy=False, verbosity=1

Project data into the frozen GDR space of reference.

Parameters

data

Query cells: an AnnData, a CytomeDataset, or a path to either. Cytome queries are scored by streaming; the expression matrix is never materialised.

reference

The reference the GDR space was built on — an AnnData/CytomeDataset/path carrying uns['gdr_reference'] / metadata['gdr_reference'] from runGDR(save_reference=True). Parts of that state are completed on first use and cached back, so the second call is free.

mode

'reference' scales the query’s score columns by the reference column norms; 'self' uses the query’s own, i.e. treats it as a new batch (what runGDR does for each reference batch when batch_key is set). Falls back to 'reference' below min_cells_self_mode cells, where self-scaling is measurably worse.

modality

Modality of the query. Cytome only — ignored for AnnData inputs, which have a single feature space (var_names).

reference_modality

Modality of the reference. Defaults to whatever the saved state recorded, so existing calls are unaffected. Supply it only to be explicit, or when the reference has no recorded modality.

These are deliberately two parameters because they are two objects. Projecting a gene-activity query onto an RNA reference is a legitimate thing to want — but doing it by accident is not, and it cannot be caught downstream: RNA and GA are both keyed by gene symbol, so the marker-recovery guard below sees matching names and passes while the values come from the wrong assay. A mismatch is therefore reported here, loudly.

layer

Layer to score. Defaults to whatever the reference used. Applies to both AnnData and cytome — there is no separate cytome_layer.

counts_layer

AnnData only. If given, the query is INFOG-normalised with the reference’s frozen constants before scoring, which is the correct choice: it keeps the query on the reference’s normalisation scale.

Notes

GDR axes are the reference’s marker-set scores, so a cell type absent from the reference still receives confident-looking coordinates — measured, a held-out subclass lands on its nearest relative (SST-Chodl→SST, PV→PV-Chandelier, VIP→SNCG 95 %) with no intrinsic warning. Read obs['<key>_novelty'] and “uns[‘_projection’][‘novelty_test’]` before interpreting the coordinates. The per-cell flag is well calibrated on false positives but low-powered (TPR 0.12 at q99); the population-level shift is the reliable readout.

queryPIASOmarkerDB

queryPIASOmarkerDB(
    gene: 'Union[str, List[str]]' = None,
    cell_type: 'Union[str, List[str]]' = None,
    study: 'str' = None,
    species: 'str' = None,
    tissue: 'str' = None,
    condition: 'str' = None,
    min_score: 'float' = None,
    max_score: 'float' = None,
    limit: 'int' = None,
    as_dict: 'Union[bool, str]' = False,
    list_studies: 'bool' = False,
    list_cell_types: 'bool' = False,
    list_genes: 'bool' = False,
)

Signature defaults

Query PIASOmarkerDB for cell type marker genes.

This is the main entry point for accessing PIASOmarkerDB through PIASO.

Parameters

gene — str or list of str, optional

Gene symbol(s) to filter by.

cell_type — str or list of str, optional

Cell type(s) to filter by.

study — str, optional

Study/publication to filter by.

species — str, optional

Species to filter by (e.g., “Human”, “Mouse”).

tissue — str, optional

Tissue to filter by.

condition — str, optional

Condition to filter by.

min_score — float, optional

Minimum specificity score (>= 0).

max_score — float, optional

Maximum specificity score (>= 0).

limit — int, optional

Maximum results to return. Default: None (no limit).

as_dict — bool or {‘both’}, optional

list_studies — bool, optional

If True, return list of available study names instead of markers. Default: False.

list_cell_types — bool, optional

If True, return list of available cell types instead of markers. Default: False.

list_genes — bool, optional

If True, return list of unique gene symbols instead of markers. Default: False.

Returns

pd.DataFrame

Marker query results (default).

dict

If as_dict=True: {cell_type: [genes]}.

tuple (pd.DataFrame, dict) If as_dict=‘both’ (deprecated).

list of str If list_studies=True, list_cell_types=True, or list_genes=True.

Examples

Query marker genes:

>>> import piaso
>>> df = piaso.tl.queryPIASOmarkerDB(gene="Foxp2", species="Mouse")
>>> df = piaso.tl.queryPIASOmarkerDB(gene=["Foxp2", "Syt6", "Tle4"])

Get the marker dictionary, ready for predictCellTypeByMarker:

>>> marker_dict = piaso.tl.queryPIASOmarkerDB(
...     study="AllenWholeMouseBrain_isocortex",
...     species="Mouse",
...     as_dict=True
... )
>>> print(f"Cell types in dict: {len(marker_dict)}")
>>> piaso.tl.predictCellTypeByMarker(adata, marker_gene_set=marker_dict)

List available studies:

>>> studies = piaso.tl.queryPIASOmarkerDB(list_studies=True)
>>> print(f"Total studies: {len(studies)}")

List cell types:

>>> cell_types = piaso.tl.queryPIASOmarkerDB(list_cell_types=True, species="Mouse")

See Also

analyzeMarkers : Analyze gene lists for cell type inference PIASOmarkerDB : Client class for advanced usage

Notes

PIASOmarkerDB website: https://piaso.org/piasomarkerdb/

read_selected_mask

_read_selected_mask(ds, modality_or_table)

Signature defaults

ds, modality_or_table

The boolean ‘highly_variable’ column of a modality’s feature table.

Takes a modality (‘ATAC’, ‘RNA’, ‘GA’, ‘tiles’) or the raw table name. Resolved through ds.features(modality) rather than getattr(ds, table): for tiles and gene activity the attribute returns a Modality, which has no __getitem__ over columns, so the old form raised or read the wrong thing for every modality but ATAC and RNA.

run_TFIDF

run_TFIDF(
    data=<unset>,
    layer: Optional[str] = None,
    scale_factor: float = 10000.0,
    streaming: bool = False,
    batch_size: int = 1024,
    output_layer: str = 'tfidf',
    measurement: Optional[str] = None,
    modality: str = 'ATAC',
    inplace: bool = False,
    source=<unset>,
    adata=<unset>,
)

Signature defaults

data=<unset>, layer=None, scale_factor=10000.0, streaming=False, batch_size=1024, output_layer='tfidf', measurement=None, modality='ATAC', inplace=False, source=<unset>, adata=<unset>

Compute TF-IDF normalization for peak count data.

Parameters

source — AnnData or CytomeDataset or str

Input data. AnnData for in-memory, Cytome for streaming.

layer — str, optional

Input layer to read counts from.

AnnData: adata.layers[layer] (or adata.X if None).
Cytome: takes precedence over measurement if both are set.

scale_factor — float, default 1e4

Scaling factor for TF values before log1p.

streaming — bool, default False

Force streaming mode even for AnnData input.

batch_size — int, default 1024

Cells per chunk in streaming mode.

output_layer — str, default ‘tfidf’

Name for the output. AnnData: writes adata.layers[output_layer]. Cytome: creates a measurement matrix named {modality}_{output_layer}. Pass output_layer=None (AnnData only) to skip the layer write and only mutate adata.X (the legacy in-place behaviour — requires inplace=True).

measurement — str, optional

Input Cytome measurement name (default: ‘counts’).

modality — str, default ‘ATAC’

Modality prefix for Cytome layer names (e.g., ‘ATAC’, ‘tiles’).

inplace — bool, default False

AnnData only. If True, also overwrite adata.X with the TF-IDF result. The default (False) writes only to adata.layers[output_layer] and leaves adata.X untouched. Useful when downstream calls (e.g. infog_svd(layer=None)) expect TF-IDF on .X.

data — AnnData, CytomeDataset or str

The data. The preferred name for adata.

Other Parameters

adata

Deprecated alias for data; passing it warns.

Returns

None

Modifies source in-place. AnnData: writes adata.layers[output_layer] (and optionally adata.X when inplace=True). Cytome: materialises {modality}_{output_layer}.

runCOSGParallel

runCOSGParallel(
    adata,
    batch_key: str,
    groupby: str = None,
    layer: str = None,
    infog_layer: str = None,
    n_svd_dims: int = 50,
    n_svd_iter: int = 7,
    n_highly_variable_genes: int = 5000,
    verbosity: int = 0,
    resolution: float = 1.0,
    mu: float = 1.0,
    n_gene: int = 30,
    use_highly_variable: bool = True,
    return_gene_names: bool = False,
    max_workers: int = 8,
    random_seed: int = 1927,
    expressed_pct: float = 0.1,
    allow_non_integer: bool = False,
)

Signature defaults

adata, batch_key, groupby=None, layer=None, infog_layer=None, n_svd_dims=50, n_svd_iter=7, n_highly_variable_genes=5000, verbosity=0, resolution=1.0, mu=1.0, n_gene=30, use_highly_variable=True, return_gene_names=False, max_workers=8, random_seed=1927, expressed_pct=0.1, allow_non_integer=False

Run COSG on batches in parallel using shared memory and multiprocessing.

Parameters

adata — AnnData

Annotated data matrix.

batch_key — str

The key in adata.obs used to identify batches.

groupby — str, optional (default: None)

The key in adata.obs used to group observations for clustering. If None, clustering will be performed.

n_svd_dims — int, optional (default: 50)

Number of SVD components to compute.

n_svd_iter — int, optional, default=7

Number of iterations for randomized SVD solver. The default is larger than the default in randomized_svd to handle sparse matrices that may have large slowly decaying spectrum. Also larger than the n_iter default value (5) in the TruncatedSVD function.

n_highly_variable_genes — int, optional (default: 5000)

Number of highly variable genes to use for SVD.

verbosity — int, optional (default: 0)

Level of verbosity for logging information.

resolution — float, optional (default: 1.0)

Resolution parameter for clustering.

layer — str, optional (default: None)

Layer of the adata object to use for COSG.

infog_layer — str, optional (default: None)

If specified, the INFOG normalization will be calculated using this layer of adata.layers, which is expected to contain the UMI count matrix. Defaults to None.

mu — float, optional (default: 1.0)

COSG parameter to control regularization.

n_gene — int, optional (default: 30)

Number of marker genes to compute for each cluster.

use_highly_variable — bool, optional (default: True)

Whether to use highly variable genes for SVD.

return_gene_names — bool, optional (default: False)

Whether to return gene names instead of indices in the marker gene DataFrame.

max_workers — int, optional (default: 8)

Maximum number of parallel workers to use. If None, defaults to the number of available CPU cores.

random_seed — int, optional

Random seed for reproducibility. Default is 1927.

Returns

DataFrame

Combined marker gene DataFrame with batch-specific suffixes.

Examples

>>> import anndata
>>> import piaso
>>> adata = anndata.read_h5ad('example_data.h5ad')
>>> marker_genes = piaso.tl.runCOSGParallel(
...     adata=adata,
...     batch_key='batch',
...     groupby=None,
...     n_svd_dims=50,
...     n_highly_variable_genes=5000,
...     verbosity=1,
...     resolution=1.0,
...     layer='log1p',
...     mu=1.0,
...     n_gene=30,
...     use_highly_variable=True,
...     return_gene_names=True,
...     max_workers=4
... )
>>> print(marker_genes.head())

runGDR

runGDR(
    data=<unset>,
    batch_key: str = None,
    groupby: str = None,
    n_gene: int = 20,
    mu: float = 10.0,
    layer: str = 'infog',
    score_layer=<default>,
    infog_layer: Optional[str] = None,
    use_highly_variable: bool = True,
    n_highly_variable_genes: int = 5000,
    n_svd_dims: int = 50,
    n_svd_iter: int = 7,
    resolution: float = 1.0,
    scoring_method: str = None,
    key_added: str = None,
    max_workers: int = 8,
    calculate_score_multiBatch: bool = True,
    n_concurrent_batches: int = None,
    verbosity: int = 0,
    random_seed: int = 1927,
    modality: str = 'RNA',
    batch_size_cytome: int = 1024,
    write_to_cytome: bool = True,
    cytome_marker_gene_key: str = 'runGDR_marker_genes',
    save_reference: bool = True,
    max_batch_cache_bytes: int = 536870912,
    expressed_pct: float = 0.1,
    allow_non_integer: bool = False,
    score_chunk_size: Optional[int] = None,
    max_score_chunk_bytes: int = 268435456,
    max_score_batch_cache_bytes: int = 536870912,
    stage1_workers: int = None,
    stage3_workers: int = None,
    adata=<unset>,
)

Signature defaults

data=<unset>, batch_key=None, groupby=None, n_gene=20, mu=10.0, layer='infog', score_layer=<default>, infog_layer=None, use_highly_variable=True, n_highly_variable_genes=5000, n_svd_dims=50, n_svd_iter=7, resolution=1.0, scoring_method=None, key_added=None, max_workers=8, calculate_score_multiBatch=True, n_concurrent_batches=None, verbosity=0, random_seed=1927, modality='RNA', batch_size_cytome=1024, write_to_cytome=True, cytome_marker_gene_key='runGDR_marker_genes', save_reference=True, max_batch_cache_bytes=536870912, expressed_pct=0.1, allow_non_integer=False, score_chunk_size=None, max_score_chunk_bytes=268435456, max_score_batch_cache_bytes=536870912, stage1_workers=None, stage3_workers=None, adata=<unset>

Run GDR (marker Gene-guided dimensionality reduction) on single-cell data.

GDR performs dimensionality reduction guided by marker genes to better preserve biological signals. When max_workers > 1 (the default), multi-batch processing uses parallel COSG marker identification and parallel gene-set scoring for faster execution. Set max_workers=1 for sequential processing (useful for debugging or memory-constrained environments).

Parameters

adata — AnnData, cytome Dataset, or str

Annotated data matrix. Also accepts a cytome Dataset or a path to a .cytome file.

batch_key — str, optional

Key in adata.obs representing batch information. Defaults to None. If provided, marker gene identifications will be performed for each batch separately.

groupby — str, optional

Key in adata.obs to specify which cell group information to use. Defaults to None. If none, de novo clustering will be performed.

n_gene — int, optional

Number of marker genes per group from COSG. Defaults to 20.

mu — float, optional

Gene expression specificity parameter, used in COSG. Defaults to 10.0.

layer — str, optional

Layer in adata.layers used for COSG marker identification. Defaults to 'infog' — PIASO’s recommended normalization for marker calling. Run piaso.tl.infog(adata) first to materialise this layer. Pass layer=None to fall back to adata.X (requires scanpy for HVG selection).

score_layer — str, optional

Layer in adata.layers used for gene-set scoring. Defaults to 'infog' (matches the recommended layer default — both COSG and score read the same INFOG-normalised matrix). Pass None to score on adata.X directly. Important: for equivalence with the cytome path, AnnData and cytome MUST score on the same data — score_layer means the same thing on both backends.

score_chunk_size — int, optional

Cytome only. Rows handed to the scoring kernel per call in stage 3. None (default) derives it from max_score_chunk_bytes, the dataset’s nonzeros per cell and the number of marker sets. Note it also blocks the first pass’s per-feature sums, so changing it perturbs the embedding in the last bits — it is not a free tuning knob.

max_score_chunk_bytes — int, default 256 MB

Cytome only. Memory budget behind score_chunk_size, TOTAL across the concurrent scoring workers (each holds one chunk). Past a few thousand rows per call the speed curve is flat, so there is little reason to raise it.

max_score_batch_cache_bytes — int, default 512 MB

Cytome only. Budget for holding a batch’s chunks between stage 3’s two streaming passes, so the second does not re-read and re-decompress rows the first just read. TOTAL across the concurrent scoring workers. A batch needs roughly n_cells * nnz_per_cell * 8 bytes; batches that do not fit stream twice, as before.

Sizing it, measured on a 200k-cell / 35-batch dataset whose batches run 384 to 13,105 cells: the 512 MB default covers 30 of the 35 and buys roughly 10 s per 100 MB until every batch fits, then nothing. Raise it if your batches are large and you have the memory; 0 disables caching.

allow_non_integer — bool, default False

INFOG refuses input whose values are not integers, because its dispersion model is defined on raw UMI counts. Set True to run on Smart-seq2 TPM/FPKM, imputed or already-corrected matrices; prefer infog_layer when raw counts do exist somewhere in the object.

infog_layer — str, optional

Source layer for piaso.tl.infog when INFOG is auto-computed (only when groupby=None and layer='infog' triggers de novo clustering). None (default) → adata.X is used as the raw-counts source. If your adata.X is normalized, point infog_layer at the layer that holds raw UMI counts (e.g. infog_layer='counts').

use_highly_variable — bool, optional

Whether to use only highly variable genes when rerunning the dimensionality reduction. Defaults to True. Only effective when groupby=None.

n_highly_variable_genes — int, optional

Number of highly variable genes to use when use_highly_variable is True. Defaults to 5000. Only effective when groupby=None.

n_svd_dims — int, optional

Number of dimensions to use for SVD. Defaults to 50. Only effective when groupby=None.

n_svd_iter — int, optional, default=7

resolution — float, optional

Resolution parameter for de novo clustering. Defaults to 1.0. Only effective when groupby=None.

scoring_method — str, optional

Specifies the gene set scoring method used to compute gene scores.

key_added — str, optional

Key under which the GDR dimensionality reduction results will be stored in adata.obsm. If None, results will be saved to adata.obsm[X_gdr].

max_workers — int, optional

Maximum number of workers for parallel computation. When > 1, multi-batch COSG and scoring run in parallel. Defaults to 8.

calculate_score_multiBatch — bool, optional

.. deprecated:: This parameter will be removed in a future version. Use max_workers=1 for sequential processing instead. Whether to calculate gene scores across multiple batches in parallel. Defaults to True.

n_concurrent_batches — int, optional

Number of batches to process concurrently. If None, auto-determined from max_workers. Default is None.

verbosity — int, optional

Verbosity level of the function. Higher values provide more detailed logs. Defaults to 0.

random_seed — int, optional

Random seed for reproducibility. Default is 1927.

modality — str, optional

Modality for cytome datasets. Defaults to 'RNA'.

layer — str, optional

Layer used for COSG marker identification, on both backends. Defaults to 'infog' (INFOG is the recommended normalization). For a cytome, run piaso.tl.infog(ds, save_layer=True) first to materialise the {modality}_infog matrix; pass layer='counts' for raw counts.

score_layer — str, optional

Layer used for gene-set scoring, on both backends. Defaults to 'infog' and mirrors layer when left unset — keep the two equal for AnnData/cytome equivalence.

batch_size_cytome — int, optional

Batch size for cytome streaming. Defaults to 1024.

write_to_cytome — bool, default True

Cytome path only. If True, the X_gdr embedding is persisted via ds.add_embedding('X_gdr', ...) and marker genes via ds.metadata[cytome_marker_gene_key]. If False, the function returns (X_gdr, marker_gene) without writing.

cytome_marker_gene_key — str, default ‘runGDR_marker_genes’

Cytome path only. Metadata key under which the marker-gene table is stored when write_to_cytome=True.

data — AnnData, CytomeDataset or str

The data. The preferred name for adata.

save_reference — bool, default True

Record the marker sets and block structure in uns['gdr_reference'] so projectGDR can map new cells into this space.

max_batch_cache_bytes — int, default 512 MB

Cytome only: budget for holding one batch’s chunks in memory between stage 3’s two passes, total across the workers.

expressed_pct — float, default 0.1

COSG’s expression fraction below which a gene is not a marker.

stage1_workers — int, optional

Cytome only: workers for the per-batch marker stage; None uses max_workers.

stage3_workers — int, optional

Cytome only: workers for the scoring stage; None uses max_workers.

Returns

None or (X_gdr, marker_gene) tuple

AnnData path: writes to adata.obsm[key_added] and returns None.
Cytome path with write_to_cytome=True (default): writes ds.embeddings['X_gdr'] + marker genes to ds.metadata and returns None.
Cytome path with write_to_cytome=False: returns (X_gdr, marker_gene) for the caller to handle.

Examples

>>> import anndata
>>> import piaso
>>>
>>> adata = anndata.read_h5ad("example.h5ad")
>>> piaso.tl.infog(adata)  # compute INFOG normalization first
>>> piaso.tl.runGDR(
...     adata,
...     batch_key="batch",
...     groupby="CellTypes",
...     n_gene=30,
...     max_workers=8,
...     verbosity=0
... )
>>> print(adata.obsm["X_gdr"])

runGDRParallel

runGDRParallel(
    data=<unset>,
    batch_key: str = None,
    groupby: str = None,
    n_gene: int = 20,
    mu: float = 10.0,
    layer: str = 'infog',
    score_layer=<default>,
    infog_layer: Optional[str] = None,
    use_highly_variable: bool = True,
    n_highly_variable_genes: int = 5000,
    n_svd_dims: int = 50,
    n_svd_iter: int = 7,
    resolution: float = 1.0,
    scoring_method: str = None,
    key_added: str = None,
    max_workers: int = 8,
    calculate_score_multiBatch: bool = True,
    n_concurrent_batches: int = None,
    verbosity: int = 0,
    random_seed: int = 1927,
    modality: str = 'RNA',
    batch_size_cytome: int = 1024,
    write_to_cytome: bool = True,
    cytome_marker_gene_key: str = 'runGDR_marker_genes',
    save_reference: bool = True,
    max_batch_cache_bytes: int = 536870912,
    expressed_pct: float = 0.1,
    allow_non_integer: bool = False,
    score_chunk_size: Optional[int] = None,
    max_score_chunk_bytes: int = 268435456,
    max_score_batch_cache_bytes: int = 536870912,
    stage1_workers: int = None,
    stage3_workers: int = None,
    adata=<unset>,
)

Signature defaults

Run GDR (marker Gene-guided dimensionality reduction) on single-cell data.

Parameters

adata — AnnData, cytome Dataset, or str

Annotated data matrix. Also accepts a cytome Dataset or a path to a .cytome file.

batch_key — str, optional

Key in adata.obs representing batch information. Defaults to None. If provided, marker gene identifications will be performed for each batch separately.

groupby — str, optional

Key in adata.obs to specify which cell group information to use. Defaults to None. If none, de novo clustering will be performed.

n_gene — int, optional

Number of marker genes per group from COSG. Defaults to 20.

mu — float, optional

Gene expression specificity parameter, used in COSG. Defaults to 10.0.

layer — str, optional

score_layer — str, optional

score_chunk_size — int, optional

max_score_chunk_bytes — int, default 256 MB

max_score_batch_cache_bytes — int, default 512 MB

allow_non_integer — bool, default False

infog_layer — str, optional

use_highly_variable — bool, optional

Whether to use only highly variable genes when rerunning the dimensionality reduction. Defaults to True. Only effective when groupby=None.

n_highly_variable_genes — int, optional

Number of highly variable genes to use when use_highly_variable is True. Defaults to 5000. Only effective when groupby=None.

n_svd_dims — int, optional

Number of dimensions to use for SVD. Defaults to 50. Only effective when groupby=None.

n_svd_iter — int, optional, default=7

resolution — float, optional

Resolution parameter for de novo clustering. Defaults to 1.0. Only effective when groupby=None.

scoring_method — str, optional

Specifies the gene set scoring method used to compute gene scores.

key_added — str, optional

Key under which the GDR dimensionality reduction results will be stored in adata.obsm. If None, results will be saved to adata.obsm[X_gdr].

max_workers — int, optional

Maximum number of workers for parallel computation. When > 1, multi-batch COSG and scoring run in parallel. Defaults to 8.

calculate_score_multiBatch — bool, optional

n_concurrent_batches — int, optional

Number of batches to process concurrently. If None, auto-determined from max_workers. Default is None.

verbosity — int, optional

Verbosity level of the function. Higher values provide more detailed logs. Defaults to 0.

random_seed — int, optional

Random seed for reproducibility. Default is 1927.

modality — str, optional

Modality for cytome datasets. Defaults to 'RNA'.

layer — str, optional

score_layer — str, optional

Layer used for gene-set scoring, on both backends. Defaults to 'infog' and mirrors layer when left unset — keep the two equal for AnnData/cytome equivalence.

batch_size_cytome — int, optional

Batch size for cytome streaming. Defaults to 1024.

write_to_cytome — bool, default True

cytome_marker_gene_key — str, default ‘runGDR_marker_genes’

Cytome path only. Metadata key under which the marker-gene table is stored when write_to_cytome=True.

data — AnnData, CytomeDataset or str

The data. The preferred name for adata.

save_reference — bool, default True

Record the marker sets and block structure in uns['gdr_reference'] so projectGDR can map new cells into this space.

max_batch_cache_bytes — int, default 512 MB

Cytome only: budget for holding one batch’s chunks in memory between stage 3’s two passes, total across the workers.

expressed_pct — float, default 0.1

COSG’s expression fraction below which a gene is not a marker.

stage1_workers — int, optional

Cytome only: workers for the per-batch marker stage; None uses max_workers.

stage3_workers — int, optional

Cytome only: workers for the scoring stage; None uses max_workers.

Returns

None or (X_gdr, marker_gene) tuple

AnnData path: writes to adata.obsm[key_added] and returns None.
Cytome path with write_to_cytome=True (default): writes ds.embeddings['X_gdr'] + marker genes to ds.metadata and returns None.
Cytome path with write_to_cytome=False: returns (X_gdr, marker_gene) for the caller to handle.

Examples

>>> import anndata
>>> import piaso
>>>
>>> adata = anndata.read_h5ad("example.h5ad")
>>> piaso.tl.infog(adata)  # compute INFOG normalization first
>>> piaso.tl.runGDR(
...     adata,
...     batch_key="batch",
...     groupby="CellTypes",
...     n_gene=30,
...     max_workers=8,
...     verbosity=0
... )
>>> print(adata.obsm["X_gdr"])

runHarmony

runHarmony(data, batch_key, use_rep='X_pca', key_added=None, random_state=0)

Signature defaults

data, batch_key, use_rep='X_pca', key_added=None, random_state=0

Run Harmony batch correction on an embedding.

Corrects batch effects in a low-dimensional embedding using the Harmony algorithm (Korsunsky et al., 2019). Uses harmonypy directly without scanpy dependency.

Parameters

data — AnnData or cytome.Dataset

Input data. For AnnData: reads embedding from obsm, stores corrected embedding in obsm. For cytome: reads/stores embeddings.

batch_key — str

Column in obs (AnnData) or cells (cytome) containing batch labels.

use_rep — str, optional (default: ‘X_pca’)

Key for the embedding to correct.

key_added — str or None, optional (default: None)

Key for the corrected embedding. If None, defaults to '{use_rep}_harmony'.

random_state — int, optional (default: 0)

Random seed for reproducibility.

Returns

np.ndarray

Corrected embedding matrix (n_cells, n_components).

runSCALAR

runSCALAR(
    adata,
    specificity_matrix: pandas.core.frame.DataFrame,
    lr_pairs: pandas.core.frame.DataFrame,
    ligand_col: str = 'ligand',
    receptor_col: str = 'receptor',
    annotation_col: Optional[str] = None,
    sender_cell_types: Optional[List[str]] = None,
    receiver_cell_types: Optional[List[str]] = None,
    n_permutations: int = 1000,
    n_nearest_neighbors: int = 100,
    layer: str = None,
    random_seed: int = 42,
    rank_by_score: bool = True,
    chunk_size: int = 50000,
    prefilter_fdr: bool = True,
    prefilter_threshold: float = 0.0,
    modality: str = 'RNA',
    expressed_pct: Optional[float] = 0.1,
    groupby: Optional[str] = None,
    statistic: str = 'product',
    enrichment_pseudocount: float = 1e-12,
)

Signature defaults

adata, specificity_matrix, lr_pairs, ligand_col='ligand', receptor_col='receptor', annotation_col=None, sender_cell_types=None, receiver_cell_types=None, n_permutations=1000, n_nearest_neighbors=100, layer=None, random_seed=42, rank_by_score=True, chunk_size=50000, prefilter_fdr=True, prefilter_threshold=0.0, modality='RNA', expressed_pct=0.1, groupby=None, statistic='product', enrichment_pseudocount=1e-12

Calculates ligand-receptor interaction scores, computes permutation-based p-values using a vectorized approach, and corrects for multiple testing using FDR for each cell type-cell type pair independently.

Args: adata: AnnData object, or a .cytome path / open cytome.Dataset. A cytome is streamed chunk by chunk to build the matched-gene null’s (mean, variance) space, so the memory cost is set by the chunk size rather than the cell count; the statistics and the resulting p-values are the same either way. modality: Cytome modality supplying expression (default 'RNA'). Ignored for AnnData. specificity_matrix: DataFrame with genes as rows, cell types as columns, and specificity scores as values. lr_pairs: DataFrame listing interacting gene pairs. ligand_col: Column name for ligands in lr_pairs. receptor_col: Column name for receptors in lr_pairs. annotation_col: Optional column in lr_pairs to carry over. sender_cell_types: List of cell types to use as senders. If None, all are used. receiver_cell_types: List of cell types to use as receivers. If None, all are used. n_permutations: Number of permutations for the null distribution. n_nearest_neighbors: Number of control genes to sample from. layer: Layer in adata to use for expression. random_seed: Seed for reproducibility. rank_by_score: If True, sorts the final output by interaction_score. chunk_size: The number of interactions to process in each vectorized chunk to manage memory usage. prefilter_fdr: If True, interactions with scores <= prefilter_threshold are excluded from FDR calculation within each group and assigned an FDR of 1.0. prefilter_threshold: The score threshold used for pre-filtering before FDR calculation. expressed_pct: Detection floor, as a fraction of a cell type’s cells in which a gene is non-zero on layer (default 0.1, CellPhoneDB’s convention). An interaction is tested only if the ligand clears the floor in the sender and the receptor in the receiver; its matched controls are drawn from genes that clear the floor in the same cell type, so the null never contains a control that is not expressed there. None disables the floor and reproduces the unfloored test exactly. Requires groupby. groupby: Column of adata.obs (or of a cytome’s cells table) holding the cell-type label of every cell. Its values must include every sender and receiver, i.e. the columns of specificity_matrix. Required when expressed_pct is set. statistic: The test statistic. 'product' (default) is spec[ligand, sender] * spec[receptor, receiver]: one very specific gene can carry a pair. 'min' is the smaller of the two, so both sides must beat the controls’ weaker sides — the “are both specific” question. Its exact null factorises into two one-sided counts and costs nothing extra. enrichment_pseudocount: Guard for log2_enrichment: the null mean is floored at this value and it is added to the observation, so a zero never reaches the log. It is on the product scale, where called nulls sit at 1e-7 to 1e-3, so the default 1e-12 changes no reported value; a larger value caps the enrichment of pairs with the sparsest nulls (1e-6 halved the median on the tutorial data).

Returns: A pandas DataFrame with interaction scores, exact p-values, FDR-corrected p-values, and two diagnostics: null_matchability (fraction of the matched sets sitting below the gene in mean space — >= 0.99 means the matched-set ceiling is active and the p-value overstates) and call_breadth (fraction of tested sender-receiver pairs in which this ligand-receptor pair is called at FDR < 0.05 — results calling most of the grid carry no group-specific information). With a detection floor the table also carries testable (both genes clear the floor in their cell types) and ligand_expressed_pct / receptor_expressed_pct (the detection fractions the floor was applied to); untestable rows have p_value = 1 and leave BH.

Three columns give the score a readable magnitude. ``ligand_specificity``
and ``receptor_specificity`` are the two inputs to the statistic.
``specificity_geomean`` is their geometric mean — a monotone function
of the product, so it ranks exactly as ``interaction_score`` does
under the default statistic, but a pair of two 0.02 specificities
reads 0.02 rather than 0.0004; size plots by it. ``null_mean`` is the
mean of the interaction's own k x k matched null, and
``log2_enrichment = log2((observed + eps) / max(null_mean, eps))`` is
how far the observation sits above what expressed,
expression-matched gene pairs achieve between these two cell types —
the one number that is comparable *across* sender–receiver pairs,
because every pair is referenced to its own null. Raw scores are
not: specificity is calibrated within a cell type, so a sharp small
subclass reaches high scores because its markers are exclusive.

Notes: Input contract. specificity_matrix must be finite and non-negative: raw COSG scores over all genes (cosg.cosg(..., n_genes_user=adata.n_vars)) are the recommended input. NaN raises — it is what cosg.indexByGene leaves for genes outside a top-N list unless set_nan_to_zero=True, and a NaN observation would otherwise score at the p-value floor. Negative values raise — COSG marks genes failing its own expressed_pct filter with −1, and two such sentinels multiply to +1, the largest product in the table. Any per-column linear rescaling (dividing a cell type’s column by a constant, e.g. its IQR) leaves every p-value and every within-pair ranking unchanged, because the observed product and its null scale together; a log or other concave transform does change them, and cosg.iqrLogNormalize (which applies log1p) is therefore not recommended as input. Scores compare within a sender–receiver pair, not across pairs.

**The detection floor is applied in the test, not in the matrix.**
Zeroing lowly-detected genes in the specificity matrix would leave
the observed pair conditioned on passing the floor while its
expression-matched controls are zeroed wherever they fall below it,
so most of the null support becomes zeros and the p-value drifts
toward "detected in more cells than its peers". Here each gene's
controls are drawn among genes detected at ``expressed_pct`` in the
same cell type, so the support stays a full ``k x k`` of real scores
and the question asked is: given both genes are expressed here, is
their specificity product unusual among expressed, expression-matched
genes?

**The p-values are exact, not sampled.** Each interaction's null is the
full ``k x k`` support of matched-ligand x matched-receptor score
products (``k = n_nearest_neighbors``), enumerated rather than drawn
from: ``p = (#{support >= observed} + 1) / (k**2 + 1)``. Because a
gene's matched set excludes the gene itself by identity, the observed
product is genuinely outside its own support, and this is the exact
p-value under exchangeability of the observed product with its matched
products — the ``+1`` is that statement, not a Monte Carlo correction.
The design follows the matched-control family of scDRS (Zhang et al.,
Nature Genetics 2022; scDRS matches by mean/variance *bins*, here
controls are the k nearest neighbours in (mean, variance) space).

**The p-value floor is ``1/(k**2+1)``** (1/10001 at the default
k=100; measured on MTG, the 100th-nearest neighbour is still matched
to a median 0.016 dex in mean — barely looser than rank 30's 0.009 —
and k=100 lowers the number of saturated matched sets).
``n_permutations`` never was and is not a way past it — only
``n_nearest_neighbors`` is. BH runs within each (sender, receiver)
group deliberately, so each pair's threshold is driven by its own
p-value distribution; the smallest attainable FDR for an isolated hit
in a group testing ``m`` interactions is ``m/(k**2+1)``, and a warning
fires when that bound exceeds 0.05.

**Pre-filtering caveat**: ``prefilter_fdr`` excludes
``interaction_score <= prefilter_threshold`` rows from BH. At the
default threshold of 0 this only removes never-rejectable tests
(a non-positive specificity product has p ≈ 1). Any threshold above 0
is selection on the test statistic itself and voids formal FDR control
among the survivors — a warning says so if you set one.

runSVD

runSVD(
    data=<unset>,
    use_highly_variable: bool = True,
    n_components: int = 50,
    random_state: Optional[int] = 10,
    scale_data: bool = False,
    n_iter=<unset>,
    key_added: str = 'X_svd',
    layer: Optional[str] = None,
    verbosity: int = 0,
    streaming: bool = False,
    batch_size: int = 1024,
    measurement: Optional[str] = None,
    oversampling: int = 10,
    modality: str = 'RNA',
    cache_chunks: bool = False,
    tfidf_params: Optional[dict] = None,
    infog_params: Optional[dict] = None,
    hvg_indices_override=None,
    selected_feature_col_name: str = 'highly_variable',
    tol: Optional[float] = None,
    method: str = 'auto',
    n_threads: int = 0,
    use_engine: bool = True,
    auto_tfidf: Optional[bool] = None,
    cell_mask=None,
    return_svd: bool = False,
    compute_on_fly: bool = True,
    source=<unset>,
    adata=<unset>,
)

Signature defaults

data=<unset>, use_highly_variable=True, n_components=50, random_state=10, scale_data=False, n_iter=<unset>, key_added='X_svd', layer=None, verbosity=0, streaming=False, batch_size=1024, measurement=None, oversampling=10, modality='RNA', cache_chunks=False, tfidf_params=None, infog_params=None, hvg_indices_override=None, selected_feature_col_name='highly_variable', tol=None, method='auto', n_threads=0, use_engine=True, auto_tfidf=None, cell_mask=None, return_svd=False, compute_on_fly=True, source=<unset>, adata=<unset>

Truncated SVD dimensionality reduction.

Supports three modes:

runSVD(adata) — standard in-memory.
runSVD(adata, streaming=True) — streaming from an in-memory AnnData.
runSVD("path.cytome") / runSVD(ds) — streaming from a cytome (self-contained: writes the embedding to the cytome and returns None).

Parameters

data — AnnData, cytome.Dataset, or str

Input. An AnnData (in-memory), an open cytome Dataset, or a path to a .cytome file (the latter two stream from disk).

use_highly_variable — bool, default True

Restrict SVD to features flagged in selected_feature_col_name.

n_components — int, default 50

Number of singular components (SVD dimensions).

random_state — int or None, default 10

Random seed for the randomized SVD solver.

scale_data — bool, default False

Z-score the (selected) features before SVD.

n_iter — int, optional

The number of rounds. For 'power' and 'krylov' without tol it is exact (7 when left unset); with tol, and for 'auto', it is the ceiling and the tolerance decides when to stop below it ('auto' left unset uses its own ceiling of 30). An explicit value is honoured by every solver.

key_added — str, default ‘X_svd’

Name for the embedding (AnnData obsm[key_added] / cytome embedding; on a cytome it is stored as {modality}_{key_added without 'X_'}).

layer — str, optional

AnnData layer to read instead of .X (in-memory path).

verbosity — int, default 0

Verbosity level.

streaming — bool, default False

Force the chunked streaming path for an in-memory AnnData. A cytome always streams; this has no effect on one.

batch_size — int, default 1024

Rows per chunk for the Python chunk reader (AnnData, or a cytome with use_engine=False). The Rust engine reads the store’s own chunks and ignores it.

measurement — str, optional

Cytome measurement/layer to read (e.g. 'counts', 'infog', 'tfidf'). Defaults to the modality’s standard layer.

oversampling — int, default 10

Extra components sampled by the randomized SVD for accuracy (solver uses n_components + oversampling).

modality — str, default ‘RNA’

Cytome modality ('RNA', 'ATAC', 'GA', 'tiles') — routes the var-entity / matrix lookups via the modality registry.

cache_chunks — bool, default False

Keep the matrix in memory after the first pass, so later passes do not read the store. The Rust engine (the default for a cytome) keeps the chunks compressed, as stored — about the size of the matrix on disk — and declines, with a warning, if that would exceed half the memory available to the process (a cluster job’s limit included). The Python chunk reader keeps them decoded, several times larger. On local disk the operating system already caches the store between passes, so this mainly helps on a network mount without the option of copying the store to local scratch.

tfidf_params — dict, optional

Apply TF-IDF inline during chunk iteration. Keys: 'cell_depth' (ndarray), 'idf' (ndarray), 'scale_factor' (float); optional 'col_mask' (bool ndarray) to slice to selected peaks. Avoids a persistent TF-IDF layer.

selected_feature_col_name — str, default ‘highly_variable’

Boolean column in the var entity (adata.var / ds.genes / ds.peaks / ds.GA_genes / ds.tiles) marking the SVD features. For ATAC/tiles, the legacy 'selected' column is auto-detected with a DeprecationWarning when the default is left untouched.

infog_params — dict, optional

The INFOG counterpart of tfidf_params: apply INFOG normalisation inline during chunk iteration instead of materialising a layer.

hvg_indices_override — ndarray, optional

Explicit column indices to run SVD on, bypassing selected_feature_col_name entirely. Used by callers that have already selected features by another route.

auto_tfidf — bool, optional

Cytome only. Apply TF-IDF on the fly: load the statistics from ds.metadata['{modality}_tfidf_params'] (or compute and cache them in one pass), with the column mask from selected_feature_col_name. Left unset it is on for modality='ATAC' or 'tiles' reading raw counts with no tfidf_params — an SVD of raw peak or tile counts is never the intended analysis — and off otherwise. Pass False to factorise raw counts anyway.

method — {‘auto’, ‘power’, ‘krylov’}, default ‘auto’

The randomized solver. 'auto' stops when the leading singular values have converged (tol=1e-3), using block Lanczos when its basis fits in memory and power iteration otherwise. 'power' is subspace iteration and 'krylov' block Lanczos regardless of memory; both run exactly n_iter rounds unless tol is set, which is the way to benchmark the effect of the round count. A round is two passes over the store, plus one pass for the sketch and one for the final projection: 2 + 2 * n_iter in all. For in-memory AnnData the solver is scikit-learn’s and this only sets the default n_iter.

tol — float, optional

Stop when the leading singular values change by less than this in relative terms, rather than always running n_iter rounds. Guarded by Aitken extrapolation of the remaining change, so a slowly converging run is not mistaken for a converged one. None (the default) runs the full n_iter. With a tolerance set, n_iter becomes the cap rather than the count.

n_threads — int, default 0

Threads for the Rust SVD pass engine. 0 uses every core this process is allowed — on a cluster, the cores the job was allocated (e.g. Slurm --cpus-per-task), not the node’s. More threads than allocated cores only makes them take turns.

use_engine — bool, default True

Use the Rust pass engine for cytome input, which reads the matrix and applies TF-IDF in one fused pass. Set False to force the Python path, which is the one to compare against when a result looks wrong.

compute_on_fly — bool, default True

Apply the normalisation inline while streaming rather than requiring a materialised layer on disk.

cell_mask — ndarray, optional

Boolean mask / sorted indices to run SVD on a cell subset (streaming / cytome paths only; returns the embedding for the masked cells).

return_svd — bool, default False

Cytome input only: by default the embedding is written to the cytome and None is returned. Pass return_svd=True to also get the in-memory (embeddings, S, Vt) tuple back.

Returns

AnnData, tuple, or None AnnData input: returns the AnnData with obsm[key_added] set (or the SVD tuple for the streaming-array path). Cytome input: writes the embedding to the cytome and returns None (self-contained); returns the (emb, S, Vt) tuple when cell_mask is set or return_svd=True.

Other Parameters

source, adata Deprecated aliases for data, kept so existing notebooks keep running; either emits a FutureWarning. They show in the signature as <unset>, which means “not supplied” — passing one along with data is an error. Use data.

Notes

On a cytome the cost is passes over the matrix, so the solver choice is the lever. Power iteration with n_iter=k costs 2k + 3 passes; method='auto' (the default) stops on a tolerance instead and typically converges in 8-12 rounds on ATAC peak matrices. The engine uses every allocated core; PIASO_SVD_PROFILE=1 in the environment prints, per pass, how many threads ran and how the time split between reading the store and computing.

Every pass re-reads the store. On a network file system that costs minutes where a node-local disk costs seconds — the same twenty-pass SVD took 5 min 13 s from a network mount and 29 s from the node’s own disk — so on a cluster copy the cytome to local scratch ($TMPDIR) first and open the copy. The engine reads without file locks, which removes the worst of it, but not the data transfer itself.

runSVDLazy

runSVDLazy(*args, **kwargs)

Signature defaults

*args, **kwargs

Deprecated alias for :func:infog_svd. Use piaso.tl.infog_svd() instead.

score

score(
    data=<unset>,
    gene_list=<unset>,
    gene_weights=None,
    n_nearest_neighbors: int = 30,
    leaf_size: int = 40,
    layer=<unset>,
    random_seed: int = 1927,
    n_ctrl_set: int = 100,
    key_added: str = None,
    compute_pvalues: bool = False,
    fallback_chunk_size: int = 10000,
    chunk_size: Optional[int] = None,
    max_workers: int = 1,
    use_rust: bool = True,
    compute_on_fly: bool = True,
    precomputed_knn: numpy.ndarray = None,
    verbosity: int = 0,
    verbose: int = None,
    modality: str = 'RNA',
    batch_size: int = 1024,
    cell_mask=None,
    score_chunk_size: Optional[int] = None,
    max_score_chunk_bytes: int = 268435456,
    max_score_batch_cache_bytes: int = 1073741824,
    pvalue_to: str = 'both',
    adata=<unset>,
)

Signature defaults

data=<unset>, gene_list=<unset>, gene_weights=None, n_nearest_neighbors=30, leaf_size=40, layer=<unset>, random_seed=1927, n_ctrl_set=100, key_added=None, compute_pvalues=False, fallback_chunk_size=10000, chunk_size=None, max_workers=1, use_rust=True, compute_on_fly=True, precomputed_knn=None, verbosity=0, verbose=None, modality='RNA', batch_size=1024, cell_mask=None, score_chunk_size=None, max_score_chunk_bytes=268435456, max_score_batch_cache_bytes=1073741824, pvalue_to='both', adata=<unset>

Compute gene-set enrichment scores for one or more gene sets, on AnnData or a Cytome dataset.

Backend / modality

Pass either an AnnData (in-memory) or a Cytome Dataset (streamed from disk). For Cytome, modality (e.g. 'RNA' / 'ATAC' / 'GA') and cytome_layer select the matrix; for AnnData, layer selects adata.layers[layer] (default 'infog').

Output

Single gene set with key_added set:

AnnData: score → adata.obs[key_added]; the full per-cell p-value table (score, pval_mc [Monte-Carlo, per cell vs its own control sets], pval [pooled empirical], plus *_FDR / nlog10_*) → adata.uns[key_added].
Cytome: score → ds.cells[key_added]; the plottable columns {key_added}_pval (Monte-Carlo), {key_added}_nlog10pval, {key_added}_FDR → ds.cells and/or the full table → ds.metadata[key_added], controlled by pvalue_to ('cells' / 'metadata' / 'both'). Always also returns (score, names, pval).

Supports two modes based on the type of gene_list:

Single gene set (list of str): Computes scores and full p-value suite (Monte Carlo, pooled empirical, FDR) for one gene set. Results are stored in adata.obs and adata.uns. Returns None.

Multiple gene sets (dict, DataFrame, or list of lists): Scores all gene sets in one batched pass using a single hstack’d sparse matmul. Optionally uses the Rust piaso_score backend for 16.8x faster matmul with 200x less RAM per thread. Returns (score_matrix, gene_set_names, pval_matrix).

Parameters

adata — AnnData

The AnnData object for the gene expression matrix.

gene_list — list of str, dict, DataFrame, or list of lists

A list of gene names (single gene set), or a dict / DataFrame / list of lists mapping gene set names to gene lists (multiple gene sets).

gene_weights — array-like or list of arrays, optional

For single gene set: a list of weights matching gene_list. For multiple gene sets: a list of weight arrays, one per gene set. If None, all genes are weighted equally. Default is None.

n_nearest_neighbors — int, optional

Number of nearest neighbors for control gene sampling. Default is 30.

leaf_size — int, optional

KDTree leaf size. Default is 40.

layer — str, optional

Layer in adata.layers to use. Default is ‘infog’.

random_seed — int, optional

Random seed for reproducibility. Default is 1927.

n_ctrl_set — int, optional

Number of control gene sets. Default is 100.

key_added — str, optional

Key for storing results in adata (single-set mode only). Default is None (‘INFOG_score’).

compute_pvalues — bool, optional

Compute Monte Carlo p-values in multi-set mode. Single-set mode always computes full p-values. Default is False.

fallback_chunk_size — int, default 10000

Rows per dense block in the pure-Python matmul fallback, which runs only when the compiled extension is unavailable. It sizes an allocation of (fallback_chunk_size, n_sets * n_ctrl_set) float64 — with 900 gene sets and n_ctrl_set=100 the default is 7.3 GB, so lower it if you are on the fallback with many gene sets. It has no effect on the Rust path, whose unit of parallel work is derived from the row count and thread count.

chunk_size — int, optional

Deprecated alias for fallback_chunk_size. It used to mean two things — this, and the Rust kernel’s tile — and the tile is now derived.

score_chunk_size — int, optional

Cytome only. Rows handed to the scoring kernel per call. None (default) derives it from max_score_chunk_bytes, the dataset’s nonzeros per cell and the number of gene sets. Note this also blocks the first pass’s per-feature sums, so changing it perturbs results in the last bits (~1e-2 on a score of order 1) — it is not a free tuning knob.

max_score_chunk_bytes — int, default 256 MB

Cytome only. Memory budget behind score_chunk_size. Bigger chunks amortise the kernel’s per-call setup; past a few thousand rows the curve is flat, so there is little reason to raise this.

max_score_batch_cache_bytes — int, default 1 GB

Cytome only. Budget for holding a batch’s chunks between the two streaming passes, so the second pass does not re-read and re-decompress rows the first pass just read (about 45% of pass 1). A batch needs roughly n_cells * nnz_per_cell * 8 bytes; batches that do not fit simply stream twice, as before.

Sizing it: on a 200k-cell dataset with 35 batches of 384-13,105 cells, the 512 MB default covers 30 of them and buys roughly 10 s per 100 MB until every batch fits, then nothing. Raise it if your batches are large and you have the memory; set it to 0 to disable caching entirely.

max_workers — int, optional

Thread count for Rust backend (1 = single-threaded). Default is 1.

use_rust — bool, optional

Try Rust fused matmul-reduce backend if available. Default is True.

precomputed_knn — ndarray, optional

Pre-computed KNN indices from _precompute_stats(). If provided, skips the KDTree construction and KNN search. Useful when calling score() multiple times on the same expression matrix.

verbosity — int, optional

Level of verbosity. Default is 0.

data — AnnData, CytomeDataset or str

The data. The preferred name for adata.

compute_on_fly — bool, default True

Cytome only: normalise the raw layer while streaming rather than requiring a materialised layer.

modality — str, default ‘RNA’

Cytome only: the modality to score.

batch_size — int, default 1024

Cytome only: rows per streaming chunk.

cell_mask — ndarray, optional

Cytome only: boolean mask or sorted indices of the cells to score; the result covers those cells only.

pvalue_to — {‘cells’, ‘metadata’, ‘both’}, default ‘both’

Cytome only: where p-values go — a cells column, the full per-cell table in ds.metadata[key_added], or both.

Other Parameters

verbose

Deprecated alias for verbosity.

Returns

Single-set mode: None. Modifies adata in-place. Multi-set mode: (score_matrix, gene_set_names, pval_matrix).

Example

>>> import piaso
>>> # Single gene set
>>> piaso.tl.score(adata, ['Gene1', 'Gene2', 'Gene3'], key_added='my_score')
>>>
>>> # Multiple gene sets (batched, with optional Rust acceleration)
>>> scores, names, pvals = piaso.tl.score(
...     adata, {'SetA': ['Gene1', 'Gene2'], 'SetB': ['Gene3', 'Gene4']},
...     compute_pvalues=True, max_workers=8
... )

Notes

The matched-control-set design and the empirical-null p-value follow scDRS (Zhang et al., Nature Genetics 54, 1572-1580, 2022); the idea is shared in spirit with AUCell and Vision. One deliberate difference in the matching itself: scDRS draws controls from mean-expression x variance bins, whereas here each gene’s controls are its k nearest neighbours in (mean, variance) space — a tighter, per-gene match. What differs here is the evaluation: all gene sets and all of their control sets are packed into one sparse weight matrix and computed as a single matrix multiplication per chunk, rather than looping and column-subsetting once per set. The loop is O(nnz) per set and grows linearly in the number of controls; the matmul is essentially flat (measured 1.95/6.45/24.41 s versus 0.55/0.62/0.78 s at 20/80/320 sets). The multiplication runs in a Rust kernel over streamed chunks, so the same path serves an in-memory AnnData and a cytome larger than RAM.

smoothCellTypePrediction

smoothCellTypePrediction(
    adata,
    groupby: str,
    use_rep: str = 'X_pca',
    k_nearest_neighbors: int = 5,
    return_confidence: bool = False,
    inplace: bool = True,
    use_existing_adjacency_graph: bool = True,
    use_faiss: bool = False,
    key_added: str = None,
    verbosity: int = 1,
    n_jobs: int = -1,
)

Signature defaults

adata, groupby, use_rep='X_pca', k_nearest_neighbors=5, return_confidence=False, inplace=True, use_existing_adjacency_graph=True, use_faiss=False, key_added=None, verbosity=1, n_jobs=-1

Smooth cell type predictions using k-nearest neighbors in a low-dimensional embedding.

Parameters

adata — AnnData

AnnData object containing single-cell data

groupby — str

Key in adata.obs containing the cell type predictions to smooth

use_rep — str, default=‘X_pca’

Key in adata.obsm containing the low-dimensional embedding to use for finding neighbors

k_nearest_neighbors — int, default=5

Number of neighbors to consider (including the cell itself)

return_confidence — bool, default=False

Whether to return confidence scores (proportion of neighbors with the majority label)

inplace — bool, default=True

Whether to modify adata inplace or return a copy

use_existing_adjacency_graph — bool, default=True

Whether to use existing neighborhood graph (adata.obsp[‘connectivities’]) if available

use_faiss — bool, default=False

Whether to use FAISS for faster neighbor search (requires faiss package)

key_added — str or None, default=None

If provided, use this key as the output key in adata.obs instead of ‘{groupby}_smoothed’

verbosity — int, default=1

Level of verbosity (0=no output, 1=basic info, 2=detailed info)

n_jobs — int, default=-1

Number of jobs for parallel processing. -1 means using all processors.

Returns

If inplace=True: None, but adds ‘groupby_smoothed’ (or key_added) to adata.obs If return_confidence=True, also adds ‘groupby_confidence’ (or key_added_confidence) to adata.obs If inplace=False: Copy of adata with added columns

Examples

>>> import piaso
>>>
>>> # Basic usage
>>> piaso.tl.smoothCellTypePrediction(
...     adata,
...     groupby='CellTypes_pred',
...     use_rep='X_pca',
...     key_added='CellTypes_pred_smoothed'
... )
>>>
>>> # With confidence scores
>>> piaso.tl.smoothCellTypePrediction(
...     adata,
...     groupby='CellTypes_pred',
...     k_nearest_neighbors=15,
...     return_confidence=True,
...     key_added='CellTypes_pred_smoothed'
... )

specificity_matrix

_specificity_matrix(
    data,
    groupby: 'str',
    modality: 'str' = 'RNA',
    cosg_mu: 'float' = 1.0,
    cosg_expressed_pct: 'float' = 0.1,
    cosg_layer: 'str' = 'counts',
    cosg_batch_size=None,
    use_cached: 'bool' = True,
    cosg_cache_key: 'str | None' = None,
    verbose: 'int' = 1,
)

Signature defaults

data, groupby, modality='RNA', cosg_mu=1.0, cosg_expressed_pct=0.1, cosg_layer='counts', cosg_batch_size=None, use_cached=True, cosg_cache_key=None, verbose=1

COSG λ specificity matrix S [features × cell_types], for AnnData or cytome. Cytome reuses the cached :func:_ensure_cosg_dense; AnnData runs cosg.cosg at full n_genes_user and pivots uns into a dense frame.

stitchSpace

stitchSpace(
    adata: anndata._core.anndata.AnnData,
    batch_key: str,
    use_rep: str = 'X_pca',
    key_added: str = 'X_stitch',
    filter_cluster_key_added: Optional[str] = None,
    filter_pruned_graph_key: Optional[str] = None,
    filter_use_global_markers: bool = False,
    filter_leiden_resolution: float = 0.5,
    filter_leiden_n_neighbors: int = 15,
    filter_n_markers: int = 50,
    filter_marker_overlap_threshold: float = 0.1,
    filter_cosg_layer: Optional[str] = None,
    filter_cosg_mu: float = 100.0,
    filter_cosg_expressed_pct: float = 0.1,
    filter_cosg_remove_lowly_expressed: bool = True,
    filter_bbknn_neighbors_within_batch: int = 3,
    filter_bbknn_trim: Optional[int] = None,
    random_state: Optional[int] = 1927,
    correction_smooth_within_batch: bool = True,
    correction_use_mutual_sqrt_weights: bool = False,
    copy: bool = False,
    verbosity: int = 0,
)

Signature defaults

adata, batch_key, use_rep='X_pca', key_added='X_stitch', filter_cluster_key_added=None, filter_pruned_graph_key=None, filter_use_global_markers=False, filter_leiden_resolution=0.5, filter_leiden_n_neighbors=15, filter_n_markers=50, filter_marker_overlap_threshold=0.1, filter_cosg_layer=None, filter_cosg_mu=100.0, filter_cosg_expressed_pct=0.1, filter_cosg_remove_lowly_expressed=True, filter_bbknn_neighbors_within_batch=3, filter_bbknn_trim=None, random_state=1927, correction_smooth_within_batch=True, correction_use_mutual_sqrt_weights=False, copy=False, verbosity=0

Performs a batch correction using a BBKNN graph that has been pruned based on marker gene overlap between batch-specific clusters. Overlap check uses local markers and optionally global markers (controlled by filter_use_global_markers).

Clusters are identified internally using Leiden and stored in adata.obs[filter_cluster_key_added]. Markers identified by COSG. Intermediate results like the compatibility ‘hypergraph’, marker gene lists, and the pruned graph structure are stored in adata.uns and adata.obsp/adata.uns.

The correction moves each cell towards the average position of its neighbors in the pruned graph.

Parameters

adata

Annotated data matrix. Needs expression data for COSG (in .X or specified layer).

batch_key

Key in adata.obs for batch information.

use_rep

Representation in adata.obsm for BBKNN, clustering, and correction (e.g., ‘X_pca’).

key_added

Base key for storing results. Corrected embedding will be in adata.obsm[key_added]. Intermediate results stored in adata.uns.

filter_cluster_key_added

Key in adata.obs where generated batch-cluster labels will be stored. If None, a default key is generated (e.g., f”{batch_key}@leiden@res{res}”).

filter_pruned_graph_key

Base key for storing the pruned graph structure in adata.obsp and adata.uns. If None, defaults to “pruned_markers”. Connectivities/distances will be stored as {filter_pruned_graph_key}_connectivities/_distances.

filter_use_global_markers

If True, run global COSG and require BOTH local AND global marker overlap for inter-batch cluster compatibility. If False (default), only local overlap is used.

filter_leiden_resolution

Resolution parameter for internal within-batch Leiden clustering.

filter_leiden_n_neighbors

KNN parameter for internal within-batch Leiden clustering’s graph.

filter_n_markers

Number of top COSG markers to compare between clusters.

filter_marker_overlap_threshold

Minimum Jaccard index for marker overlap to consider clusters compatible.

filter_cosg_layer

Layer in adata.layers to use for COSG marker identification. If None (default), uses adata.X.

filter_cosg_mu

mu parameter for COSG (default: 100.0). Higher values increase sparsity.

filter_cosg_expressed_pct

expressed_pct parameter for COSG (default: 0.1). Minimum expression pct for a gene.

filter_cosg_remove_lowly_expressed

remove_lowly_expressed parameter for COSG (default: True). Filter lowly expressed genes.

filter_bbknn_neighbors_within_batch

neighbors_within_batch parameter for the initial bbknn.bbknn call.

filter_bbknn_trim

Optional trim parameter passed to the initial bbknn.bbknn call.

random_state

Seed for the random number generator used in Leiden clustering for reproducibility. Default: 1927.

correction_smooth_within_batch

If True, smooth the correction vector within batches using the pruned graph structure.

correction_use_mutual_sqrt_weights

If True, applies symmetrization and sqrt weighting to the pruned graph before the correction step.

copy

If True, return a modified copy of adata. Otherwise, modify adata inplace.

verbosity

Level of detail to print: 0 (minimal), 1 or higher (more progress messages and intermediate storage locations). Default: 0. Controls BBKNN logging level.

Returns

AnnData or None If copy=True, returns the modified AnnData object. Otherwise, modifies the input adata object inplace and returns None. Adds/updates:

adata.obsm[key_added]: The corrected embedding.
adata.obs[filter_cluster_key_added]: Generated batch-cluster labels (using ’@’ delimiter).
adata.uns[f'{key_added}_hypergraph_compatibility']: Compatibility dict.
adata.uns[f'{key_added}_local_markers']: Local marker dict.
adata.uns[f'{key_added}_global_markers']: Global marker dict.
adata.obsp[f'{filter_pruned_graph_key}_connectivities']: Pruned graph connectivities.
adata.obsp[f'{filter_pruned_graph_key}_distances']: Pruned graph dummy distances.
adata.uns[filter_pruned_graph_key]: Neighbors dictionary for pruned graph.
adata.uns[f'{key_added}_params']: Dictionary of parameters used.

Example

>>> import anndata
>>> import piaso
>>> adata = anndata.read_h5ad('pbmc68k.h5ad')
>>> # Simulate batches (replace with actual batch info)
>>> adata.obs['batch'] = ['A' if i % 2 == 0 else 'B' for i in range(adata.n_obs)]
>>> # Assume normalized data is in adata.layers['log1p']
>>> adata.layers['log1p'] = adata.X.copy()
>>> # Run correction using log1p layer for COSG, increased verbosity
>>> piaso.tl.stitchSpace(
...     adata,
...     batch_key='batch',
...     use_rep='X_pca',
...     key_added='X_stitch_corrected',
...     filter_cluster_key_added='batch@cluster_stitch',
...     filter_cosg_layer='log1p',
...     random_state=1927,
...     verbosity=1
... )
>>> # Visualize results
>>> piaso.tl.neighbors(adata, use_rep='X_stitch_corrected')
>>> piaso.tl.umap(adata)
>>> piaso.pl.plotEmbedding(adata, color='batch')
>>> piaso.pl.plotEmbedding(adata, color='batch@cluster_stitch')

umap

umap(
    data,
    use_rep=None,
    min_dist=0.5,
    spread=1.0,
    n_components=2,
    random_state=42,
    key_added='X_umap',
    knn_result=None,
    neighbors_key='neighbors',
    modality='RNA',
)

Signature defaults

data, use_rep=None, min_dist=0.5, spread=1.0, n_components=2, random_state=42, key_added='X_umap', knn_result=None, neighbors_key='neighbors', modality='RNA'

Compute UMAP embedding from precomputed kNN graph.

Requires piaso.tl.neighbors() to have been run first. Pass the dict returned by neighbors() as knn_result for cytome mode.

Parameters

data — AnnData or cytome.Dataset

If AnnData: reads kNN from uns, stores UMAP in obsm. If cytome.Dataset: reads embedding, stores UMAP embedding in cytome.

use_rep — str

Name of the representation whose ARRAY is fed to umap-learn. What it controls depends on whether a stored kNN graph is being reused:

cytome, or AnnData with a stored piaso kNN: neighbors are not re-run — the embedding’s topology comes wholly from the stored graph, and use_rep only selects the data array (umap-learn’s API requires one even with precomputed_knn). Naming a different representation than the graph was built on raises, because it would change nothing while implying it had.
AnnData with no stored kNN: umap-learn computes its own neighbors from this array (its defaults, not ours) — announced with a warning. Run :func:piaso.tl.neighbors first to control and reuse the graph.

min_dist — float

Minimum distance parameter for UMAP.

spread — float

Spread parameter for UMAP.

n_components — int

Number of UMAP dimensions.

random_state — int

Random seed for reproducibility.

key_added — str

Key/name for the UMAP coordinates.

knn_result — dict, optional

Result dict from neighbors() with ‘knn_indices’ and ‘knn_dists’. Only needed for the in-memory ndarray / cell_mask path. For a cytome.Dataset it is not required — the kNN is reconstructed from the persisted distances graph (self-contained). For AnnData it falls back to uns.

neighbors_key — str

Which stored neighbors graph to reuse — matches the key_added passed to :func:piaso.tl.neighbors. Default 'neighbors' (the un-prefixed graph). Honoured on both the AnnData and cytome paths; naming a graph that was never computed raises rather than quietly falling back to the default one.

modality — str, default ‘RNA’

Cytome only: the modality the graph was stored under, and the one the embedding is written under.

Returns

None, or np.ndarray for an in-memory call Nothing, on both backends: the embedding is written to adata.obsm[key_added] or to the cytome under key_added, and read back from there. Only the in-memory form — data=None with a knn_result — returns the coordinates, because there is nowhere to write them.

Changed in 1.2.4: the AnnData path used to return the coordinates as well as writing them, so the same call had two different contracts depending on the backend and echoed a 2-column array in a notebook. Read them from adata.obsm[key_added].

Moved to cytorete

These names still work, but the method they call now lives in cytorete — pip install cytorete, then use it directly as cytorete.tl.<name>. Calling them through PIASO raises a pointer to that package if it is not installed.

inferGRN, inferRegulon, inferTFActivity, regulonActivity, regulonSpecificity