diff --git a/CHANGELOG.md b/CHANGELOG.md index c2e7f2232..4d7fb2d04 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -8,6 +8,7 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0 ## [Unreleased] ### Added +- `single-cell-pna demux`, `collapse`, and `graph` accept repeated `--panel` options. The panels are concatenated, and `marker_id`, `sequence_1`, and `sequence_2` must be unique across the result. - `pixelator.pna.analysis.differential_abundance` to run marker differential abundance (Wilcoxon via `scanpy.tl.rank_genes_groups`). Effect size is a mean difference, not scanpy's log-fold change. - Samples with no cells above the size threshold no longer fail the pipeline. A null pxl file with a reason is written and passed through later steps, so the sample still shows up downstream. @@ -26,6 +27,8 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0 - `pixelator.pna.analysis.filter_proximity_scores` to filter a proximity score table by marker abundance (Python analog of pixelatorR `FilterProximityScores`). ### Changed +- New `.pxl` files store the antibody panel in the `panels` and `panel_sources` tables. Rewriting a pixelator 0.22.0–0.30.0 file stores that panel there as well. Panel columns are joined onto the AnnData in memory when the file is read, and are not written to `var`. `uns['panel_metadata']` is no longer written, and pixel-file metadata no longer includes `panel_name` or `panel_version`. +- When several `.pxl` files are read together, each source is updated to the newest patch those files carry, matching on name, product, and major.minor. The in-memory edgelist, proximity, layouts, and `var` use the updated marker ids. The files on disk are not rewritten. - Updated pixelgen-pixelator-core to 0.2.0 improving peak memory usage in the graph step by ~20%. - `density_scatter_plot` now lives in `pixelator.plot` (previously `pixelator.mpx.plot`). - `uei_count` is now optional on PNA edgelists in `sample_calling` and the graph component @@ -38,6 +41,7 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0 when one exists. ### Removed +- `PNAAntibodyPanel.from_adata`. Panels are read from the `.pxl` file. Files from pixelator 0.22.0 through 0.30.0 are still read from `uns['panel_metadata']` and `var`. A file from before 0.22.0 has no panel; denoise and sample calling ask for the sample to be rerun with a current version of the software. - Stale nf-core/pixelator integration tests (`task test-nf-core-pixelator` and `.github/workflows/nf-core-pixelator-tests.yml`). - The `fruchterman_reingold`, `fruchterman_reingold_3d`, `kamada_kawai`, and @@ -54,6 +58,8 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0 passed directly to native community detection, and the recovered components are unchanged. ### Fixed +- `Edgelist.to_record_batches` renames legacy `marker1` and `marker2` columns to `marker_1` and `marker_2`, matching `to_polars`, `to_df`, and `iterator`. The stream on dev left those column names unchanged. +- Denoise rebuilds the count matrix from the markers already in `var`, so hashing clones removed by sample calling do not come back. The full panel is still stored in the panel tables. - `coarsened_pmds_layout` sizes PMDS pivots from the full graph when Leiden yields too few communities, so a valid low `pivots` no longer fails `pmds_layout`'s `0.2 * n` lower bound. diff --git a/src/pixelator/pna/analysis/denoise.py b/src/pixelator/pna/analysis/denoise.py index 67c13fe46..305d8e2dc 100644 --- a/src/pixelator/pna/analysis/denoise.py +++ b/src/pixelator/pna/analysis/denoise.py @@ -19,8 +19,6 @@ from pixelator.pna.analysis_engine import PerComponentTask from pixelator.pna.anndata import add_missing_adata_info, pna_edgelist_to_anndata -from pixelator.pna.config import pna_config -from pixelator.pna.config.panel import PNAAntibodyPanel, load_antibody_panel from pixelator.pna.graph import PNAGraph from pixelator.pna.graph.adaptive_core_expansion import adaptive_core_expansion from pixelator.pna.graph.node_pls import ( @@ -29,7 +27,7 @@ node_pls, ) from pixelator.pna.pixeldataset import PNAPixelDataset, read -from pixelator.pna.pixeldataset.io import PixelFileWriter, PxlFile +from pixelator.pna.pixeldataset.io import PixelFileWriter, PxlFile, read_dataset_panel logger = logging.getLogger(__name__) @@ -644,15 +642,7 @@ def add_to_pixel_file(self, data: pd.DataFrame, pxl_file_target: PxlFile): """Filter edgelist by removed nodes and write denoise metrics to AnnData.""" pxl = PNAPixelDataset.from_files(pxl_file_target) old_adata = pxl.adata() - try: - panel = PNAAntibodyPanel.from_pxl_dataset(read(pxl_file_target.path)) - except KeyError: - # If pxl file does not contain panel data, try to load it from - # the panel name. - # This will happen when old pxl files generated before v0.22.0 - # are used. - panel_name = pxl.metadata().popitem()[1]["panel_name"] - panel = load_antibody_panel(pna_config, panel_name) + panel = read_dataset_panel(read(pxl_file_target.path)) nodes_to_remove = ( data.loc[~data["umi"].isna(), "umi"].astype(np.uint64).tolist() ) @@ -662,7 +652,14 @@ def add_to_pixel_file(self, data: pd.DataFrame, pxl_file_target: PxlFile): write_denoised_edgelist(pxl, nodes_to_remove, denoised_edgelist_path) with PixelFileWriter(pxl_file_target.path) as writer: writer.write_edgelist(Path(denoised_edgelist_path)) - adata = pna_edgelist_to_anndata(writer.get_connection(), panel) + # Rebuild counts from the markers already in var. Using the full + # panel would put hashing clones removed by sample calling back. + adata = pna_edgelist_to_anndata( + writer.get_connection(), + panel, + markers=list(old_adata.var_names), + ) + writer.write_panel(panel) old_adata.obs.rename( columns={"isotype_fraction": "pre_denoise_isotype_fraction"}, inplace=True, diff --git a/src/pixelator/pna/anndata.py b/src/pixelator/pna/anndata.py index 94815b9f7..9cb09a71e 100644 --- a/src/pixelator/pna/anndata.py +++ b/src/pixelator/pna/anndata.py @@ -30,18 +30,10 @@ def calculate_antibody_metrics(counts_df): return pd.concat([total_antibody, relative_antibody, components_detected], axis=1) -def add_panel_information(adata: AnnData, panel: PNAAntibodyPanel) -> AnnData: - """Add panel data to var.""" - adata.var = adata.var.join(panel.df, how="left") - - adata.uns["panel_metadata"] = panel.metadata.model_dump() - adata.uns["panel_metadata"]["panel_columns"] = list(panel.df.columns) - - return adata - - def pna_edgelist_to_anndata( - pixel_connection: duckdb.DuckDBPyConnection, panel: PNAAntibodyPanel + pixel_connection: duckdb.DuckDBPyConnection, + panel: PNAAntibodyPanel, + markers: list[str] | None = None, ) -> AnnData: """Build an AnnData object from a DuckDB connection to a pixel file and a panel object. @@ -49,6 +41,9 @@ def pna_edgelist_to_anndata( pixel_connection: A DuckDB connection to a pixel file. The connection must contain an 'edgelist' table with the required columns (e.g., component, marker_1, marker_2, umi1, umi2, read_count). panel: The antibody panel object containing marker metadata. + markers: Marker ids to use as the count-matrix columns. Defaults to + every marker on ``panel``. Denoise passes the markers already in + ``var`` so hashing clones removed by sample calling are not added back. Returns: An AnnData object with counts and panel information. @@ -73,10 +68,11 @@ def pna_edgelist_to_anndata( .tolist() ) + marker_list = list(panel.markers if markers is None else markers) n_components = len(components) - n_markers = len(panel.markers) + n_markers = len(marker_list) component_to_idx = {c: i for i, c in enumerate(components)} - marker_to_idx = {m: i for i, m in enumerate(panel.markers)} + marker_to_idx = {m: i for i, m in enumerate(marker_list)} X = np.zeros((n_components, n_markers), dtype=np.uint32) n_umi1_arr = np.zeros(n_components, dtype=np.uint64) @@ -152,7 +148,7 @@ def pna_edgelist_to_anndata( node_counts_df = pd.DataFrame( X, index=component_index, - columns=pd.Index(panel.markers, name="marker_id"), + columns=pd.Index(marker_list, name="marker_id"), ) logger.debug("Computing component metrics.") @@ -182,7 +178,7 @@ def pna_edgelist_to_anndata( logger.debug("Computing antibody metrics.") antibody_metrics_df = calculate_antibody_metrics(counts_df=node_counts_df) - antibody_metrics_df = antibody_metrics_df.reindex(index=panel.markers, fill_value=0) + antibody_metrics_df = antibody_metrics_df.reindex(index=marker_list, fill_value=0) antibody_metrics_df.index.name = "marker_id" # Do a dtype conversion of the columns here since AnnData cannot handle # a pyarrow arrays. @@ -198,10 +194,17 @@ def pna_edgelist_to_anndata( var=antibody_metrics_df, ) - adata = add_panel_information(adata, panel) - total_marker_counts = node_counts_df.sum(axis=1) - isotype_markers = adata.var[adata.var["control"]].index + control = panel.df["control"] + if pd.api.types.is_bool_dtype(control): + control_mask = control.fillna(False).astype(bool) + else: + control_mask = control.astype(str).str.lower().eq("yes") + isotype_markers = [ + marker + for marker in control_mask.index[control_mask] + if marker in node_counts_df.columns + ] isotype_counts = node_counts_df[isotype_markers].sum(axis=1) adata.obs["isotype_fraction"] = isotype_counts / total_marker_counts diff --git a/src/pixelator/pna/cli/collapse.py b/src/pixelator/pna/cli/collapse.py index 13f4887ed..8e1475594 100644 --- a/src/pixelator/pna/cli/collapse.py +++ b/src/pixelator/pna/cli/collapse.py @@ -29,7 +29,8 @@ ) from pixelator.pna.collapse.paired.collapser import MoleculeCollapser from pixelator.pna.collapse.utilities import split_collapse_inputs -from pixelator.pna.config import load_antibody_panel, pna_config +from pixelator.pna.config import pna_config +from pixelator.pna.config.panel import load_antibody_panels from pixelator.pna.utils import get_demux_filename_info logger = logging.getLogger("collapse") @@ -124,7 +125,7 @@ def collapse( sanity_check_inputs(input_files=input_files, allowed_extensions=("parquet",)) assay = pna_config.get_assay(design) - panel = load_antibody_panel(pna_config, panel) + panel = load_antibody_panels(pna_config, panel) # create the output directory collapse_output = create_output_stage_dir(output, "collapse") diff --git a/src/pixelator/pna/cli/common.py b/src/pixelator/pna/cli/common.py index 4b06b3cee..e03885076 100644 --- a/src/pixelator/pna/cli/common.py +++ b/src/pixelator/pna/cli/common.py @@ -5,6 +5,7 @@ import functools import logging +import os from pathlib import Path import click @@ -103,16 +104,7 @@ def wrapper(*args, **kwargs): return wrapper -def validate_panel(ctx, param, value): - """Validate the panel commandline option. - - Args: - ctx: The click context - param: The click parameter - value: The click value - Returns: - The validated value - """ +def _validate_one_panel(value: str) -> str: try: if Path(value).exists(): return value @@ -133,6 +125,26 @@ def validate_panel(ctx, param, value): return value +def validate_panel(ctx, param, value): + """Validate the panel commandline option. + + ``--panel`` may be repeated. Each value is a supported panel name or a path + to a panel CSV. + + Args: + ctx: The click context + param: The click parameter + value: The click value + Returns: + The validated value + """ + if value is None: + return value + if isinstance(value, (str, os.PathLike)): + return _validate_one_panel(value) + return tuple(_validate_one_panel(item) for item in value) + + def panel_option(func): """Decorate a click command and add the --panel option.""" from pixelator.pna.config import pna_config @@ -144,10 +156,13 @@ def panel_option(func): @click.option( "--panel", required=True, - default=None, + multiple=True, type=click.UNPROCESSED, callback=validate_panel, - help="The name of a panel to load from the supported panels. Optionally, provide a path to a custom panel file.", + help=( + "The name of a panel to load from the supported panels, or a path to a " + "custom panel file. Repeat to concatenate several panels." + ), ) @functools.wraps(func) def wrapper(*args, **kwargs): diff --git a/src/pixelator/pna/cli/demux.py b/src/pixelator/pna/cli/demux.py index 3d1efa53d..773cffee4 100644 --- a/src/pixelator/pna/cli/demux.py +++ b/src/pixelator/pna/cli/demux.py @@ -25,7 +25,7 @@ threads_option, ) from pixelator.pna.config import pna_config -from pixelator.pna.config.panel import load_antibody_panel +from pixelator.pna.config.panel import load_antibody_panels from pixelator.pna.demux import ( correct_marker_barcodes, demux_barcode_groups, @@ -143,7 +143,7 @@ def demux( # load assay design assay = pna_config.get_assay(design) # load marker panel - panel = load_antibody_panel(pna_config, panel) + panel = load_antibody_panels(pna_config, panel) logger.info(f"Correcting marker barcodes for input: {fastq_file}") diff --git a/src/pixelator/pna/cli/graph.py b/src/pixelator/pna/cli/graph.py index d22bcf416..3c3ec6374 100644 --- a/src/pixelator/pna/cli/graph.py +++ b/src/pixelator/pna/cli/graph.py @@ -18,7 +18,8 @@ write_parameters_file, ) from pixelator.pna.cli.common import logger, output_option, panel_option -from pixelator.pna.config import load_antibody_panel, pna_config +from pixelator.pna.config import pna_config +from pixelator.pna.config.panel import load_antibody_panels from pixelator.pna.graph.community_detection import ( RefinementOptions, StagedRefinementOptions, @@ -257,7 +258,7 @@ def graph( ) output_path = graph_output / f"{sample_name}.graph.pxl" - panel = load_antibody_panel(pna_config, panel) + panel = load_antibody_panels(pna_config, panel) initial_stage_refinement_options = RefinementOptions( leiden_resolution=initial_stage_leiden_resolution, max_edges_to_remove=initial_stage_max_edges_to_remove, @@ -316,8 +317,6 @@ def graph( output_path, sample_name=sample_name, reason=str(exc), - panel_name=panel.name, - panel_version=panel.version, ) report = _graph_report( sample_name, diff --git a/src/pixelator/pna/cli/sample_calling.py b/src/pixelator/pna/cli/sample_calling.py index df24c2c7f..204bc5f7a 100644 --- a/src/pixelator/pna/cli/sample_calling.py +++ b/src/pixelator/pna/cli/sample_calling.py @@ -19,9 +19,8 @@ ) from pixelator.pna import read from pixelator.pna.cli.common import output_option -from pixelator.pna.config.panel import PNAAntibodyPanel from pixelator.pna.pixeldataset import NullPxlFileError -from pixelator.pna.pixeldataset.io import PxlFile, write_null_pxl +from pixelator.pna.pixeldataset.io import PxlFile, read_dataset_panel, write_null_pxl from pixelator.pna.sample_calling import ( create_final_report, sample_calling, @@ -106,7 +105,7 @@ def sample_calling_cli( undetermined_sample_name = f"{pool_name}_undetermined" try: - panel_info = PNAAntibodyPanel.from_pxl_dataset(read(input_pxl_file)) + panel = read_dataset_panel(read(input_pxl_file)) except NullPxlFileError as exc: logger.warning("%s", exc) _pass_through_null_sample_calling( @@ -118,9 +117,13 @@ def sample_calling_cli( pool_metadata=PxlFile(Path(input_pxl_file)).metadata(), ) return - hashing_antibodies_in_panel = set( - panel_info.df[panel_info.df["sample_hashing"] == "yes"].index.to_list() - ) + if "sample_hashing" not in panel.df.columns: + raise ValueError( + "Sample calling requires a sample_hashing column on the panel " + "so hashing markers can be identified. This panel has no " + "sample_hashing column." + ) + hashing_antibodies_in_panel = panel.hashing_marker_ids samplesheet_df = pl.read_csv(samplesheet) _reject_reserved_samplesheet_names( samplesheet_df["sample"].to_list(), diff --git a/src/pixelator/pna/config/__init__.py b/src/pixelator/pna/config/__init__.py index 17225b041..4dbd8afce 100644 --- a/src/pixelator/pna/config/__init__.py +++ b/src/pixelator/pna/config/__init__.py @@ -19,6 +19,7 @@ from pixelator.pna.config.panel import ( PNAAntibodyPanel, load_antibody_panel, + load_antibody_panels, ) __all__ = [ @@ -35,4 +36,5 @@ "PNAAntibodyPanel", "AntibodyPanelMetadata", "load_antibody_panel", + "load_antibody_panels", ] diff --git a/src/pixelator/pna/config/config_class.py b/src/pixelator/pna/config/config_class.py index a940e0581..55dd949bc 100644 --- a/src/pixelator/pna/config/config_class.py +++ b/src/pixelator/pna/config/config_class.py @@ -93,16 +93,23 @@ def load_panel_file(self, path: PathType) -> None: def add_panel(self, panel: PNAAntibodyPanel) -> None: """Register a panel and its lookup keys in the config. - The panel is indexed by panel name (or filename fallback), optional product, - and aliases. + The panel is indexed by panel name, optional product, and aliases. + Only a single-source panel can be registered. Panel CSV files are + one source. Args: panel: Panel object to add. Raises: + ValueError: If ``panel`` has more or fewer than one source. PanelException: If an alias already maps to a different panel key. """ - key = panel.name if panel.name is not None else str(panel.filename) + if len(panel.sources) != 1: + raise ValueError( + "Only a single-source panel can be added to the config. " + f"Got {len(panel.sources)} sources." + ) + key = panel.name self.panels[key].append(panel) # allow to also get panel by product name if provided in the panel file @@ -279,18 +286,17 @@ def get_panel( + "alias to disambiguate.", ) - def keyfunc(p): + def keyfunc(p: PNAAntibodyPanel): """Keyfunc. Args: - p: p. + p: PNAAntibodyPanel. + + Returns: + Parsed semantic version of the panel. """ version = p.version - if version is None: - v = semver.Version.parse("0.0.0") - - v = semver.Version.parse(version) - return v + return semver.Version.parse(version) panels_with_key = sorted(panels_with_key, key=keyfunc, reverse=True) return panels_with_key[0] diff --git a/src/pixelator/pna/config/panel.py b/src/pixelator/pna/config/panel.py index cf89c7a1e..eb6c2d873 100644 --- a/src/pixelator/pna/config/panel.py +++ b/src/pixelator/pna/config/panel.py @@ -5,12 +5,14 @@ from __future__ import annotations +import json +import os import warnings +from collections import defaultdict +from dataclasses import dataclass from functools import cached_property from pathlib import Path -from typing import TYPE_CHECKING, List, Optional, Set - -from anndata import AnnData +from typing import TYPE_CHECKING, List, Optional, Sequence, Set try: from typing import Self @@ -21,6 +23,8 @@ import pandas as pd import polars as pl +from anndata import AnnData +from packaging.version import Version from pixelator.common.config.panel import ( AntibodyPanelMetadata, @@ -31,7 +35,52 @@ if TYPE_CHECKING: from pixelator.pna.config.config_class import PNAConfig - from pixelator.pna.pixeldataset.dataset import PNAPixelDataset + +# Trailing ``-`` is the hash group (``B2M-1`` → ``B2M``). The same +# pattern matches ordinary names such as ``PD-1``, so it is only applied to +# rows already flagged by ``sample_hashing``. +_HASHING_MARKER_ID_RE = re.compile(r"^(?P.+)-(?P\d+)$") + + +@dataclass(frozen=True) +class PanelSource: + """One panel file that contributed markers to a ``PNAAntibodyPanel``.""" + + metadata: AntibodyPanelMetadata + file_name: str | None = None + filepath: str | None = None + + +def sample_hashing_mask(sample_hashing: pd.Series) -> pd.Series: + """Return a boolean mask for values that flag a hashing marker.""" + if pd.api.types.is_bool_dtype(sample_hashing): + return sample_hashing.fillna(False).astype(bool) + if pd.api.types.is_numeric_dtype(sample_hashing): + return sample_hashing.fillna(0).astype(bool) + normalized = sample_hashing.astype(str).str.strip().str.lower() + return normalized.isin(["yes", "true"]) + + +def split_hashing_marker_id(marker_id: str) -> tuple[str, str] | None: + """Return ``(base, index)`` for a hashing id such as ``B2M-1``.""" + match = _HASHING_MARKER_ID_RE.fullmatch(str(marker_id)) + if match is None: + return None + return match.group("base"), match.group("index") + + +def collapsed_hashing_marker_id(marker_id: str) -> str: + """Return the marker id sample calling stores for a hashing antibody.""" + parts = split_hashing_marker_id(marker_id) + return parts[0] if parts is not None else str(marker_id) + + +def _hashing_marker_ids(panel_df: pd.DataFrame) -> set[str]: + """Return hashing marker ids, or an empty set when the column is absent.""" + if "sample_hashing" not in panel_df.columns: + return set() + mask = sample_hashing_mask(panel_df["sample_hashing"]) + return {str(marker_id) for marker_id in panel_df.index[mask]} class PNAAntibodyPanel: @@ -52,27 +101,60 @@ class PNAAntibodyPanel: def __init__( self, df: pd.DataFrame, - metadata: AntibodyPanelMetadata, + metadata: AntibodyPanelMetadata | None = None, file_name: Optional[str] = None, filepath: Optional[PathType] = None, + *, + sources: list[PanelSource] | None = None, + marker_source_ids: pd.Series | None = None, ) -> None: - """Load a panel from a dataframe and metadata. + """Build a panel from a marker table. + + Pass ``metadata`` or ``sources``, not both. ``metadata`` becomes one + source, using ``file_name`` and ``filepath``. Pass ``sources`` for a + panel that already names its files. The ``metadata`` property returns + that entry only when the panel has one source. Args: - df: The dataframe containing the panel information. - metadata: The metadata for the panel. - file_name: The optional basename of the file from which the panel is loaded. - filepath: The optional full path of the file from which the panel is loaded. + df: Marker table, indexed by marker id. + metadata: Metadata for the single source. Omit when passing + ``sources``. + file_name: Basename of the file this panel was loaded from. + filepath: Full path of the file this panel was loaded from. + sources: Panel files that contributed markers. Omit when passing + ``metadata``. + marker_source_ids: Source index for each marker. Required when + there is more than one source. - Returns: - None Raises: - AssertionError: exception if panel file is missing, invalid or with incorrect format + ValueError: If both ``metadata`` and ``sources`` are omitted or + both are given, or if several sources are given without + ``marker_source_ids``. + AssertionError: If the marker table fails panel validation. """ self._filename = file_name self._filepath: Optional[Path] = Path(filepath).resolve() if filepath else None - self.metadata = metadata self._df = df + if sources is not None and metadata is not None: + raise ValueError("Pass metadata or sources, not both.") + if sources is None: + if metadata is None: + raise ValueError("Pass metadata or sources.") + sources = [ + PanelSource( + metadata=metadata, + file_name=file_name, + filepath=str(self._filepath) if self._filepath else None, + ) + ] + self.sources: list[PanelSource] = list(sources) + if marker_source_ids is None: + if len(self.sources) > 1: + raise ValueError( + "marker_source_ids is required when a panel has multiple sources." + ) + marker_source_ids = pd.Series(0, index=df.index, dtype="int64") + self._marker_source_ids = marker_source_ids # validate the panel errors = self.validate_antibody_panel(df) @@ -112,103 +194,260 @@ def from_csv(cls, filename: PathType) -> Self: return cls(df, metadata, file_name=panel_file.name, filepath=panel_file) @classmethod - def from_pxl_dataset( - cls, pxl_data: PNAPixelDataset, file_name: Optional[str] = None + def from_legacy_var( + cls, + var: pd.DataFrame, + panel_metadata: dict, + *, + file_name: str | None = None, ) -> Self: - """Create an AntibodyPanel from a pxl dataset. + """Build a panel from a pixelator 0.22.0 through 0.30.0 ``var`` table. - Args: - pxl_data: A PNAPixelDataset object. - file_name: The optional name of the file from which the pxl dataset was loaded. + ``panel_metadata`` is the ``uns['panel_metadata']`` entry, and ``var`` + is indexed by marker id. Those releases wrote ``panel_columns`` and + the columns it names together. + """ + df = var[list(panel_metadata["panel_columns"])] + metadata = AntibodyPanelMetadata.model_validate(panel_metadata) + return cls(df, metadata, file_name=file_name) - Returns: - The AntibodyPanel object. (AntibodyPanel) + @classmethod + def _from_panel_tables(cls, markers: pd.DataFrame, sources: pd.DataFrame) -> Self: + """Build a panel from the stored marker and source tables.""" + panel_sources: list[PanelSource] = [] + for row in sources.sort_values("source_id").itertuples(index=False): + aliases = json.loads(row.aliases) if isinstance(row.aliases, str) else [] + metadata = AntibodyPanelMetadata( + name=row.name, + version=row.version, + product=None if pd.isna(row.product) else row.product, + description=None if pd.isna(row.description) else row.description, + aliases=aliases or [], + archived=bool(row.archived) if not pd.isna(row.archived) else False, + ) + file_name = None if pd.isna(row.file_name) else row.file_name + filepath = None if pd.isna(row.filepath) else row.filepath + panel_sources.append( + PanelSource(metadata=metadata, file_name=file_name, filepath=filepath) + ) - Raises: - KeyError: exception if panel information is missing in the pxl dataset, - """ - logger.debug("Creating Antibody panel from PNAPixelDataset object") - adata = pxl_data.adata() - panel = cls.from_adata(adata, file_name=file_name) - logger.debug("Antibody panel from PNAPixelDataset created") - return panel + markers = markers.sort_values("row_nr") + source_ids = markers["source_id"].astype(int) + drop_cols = ["row_nr", "source_id"] + df = markers.drop(columns=[col for col in drop_cols if col in markers.columns]) + df = df.set_index("marker_id") + df.index.name = "marker_id" + source_ids.index = df.index + file_name = panel_sources[0].file_name if len(panel_sources) == 1 else None + filepath = panel_sources[0].filepath if len(panel_sources) == 1 else None + return cls( + df, + file_name=file_name, + filepath=filepath, + sources=panel_sources, + marker_source_ids=source_ids.astype("int64"), + ) @classmethod - def from_adata(cls, adata: AnnData, file_name: Optional[str] = None) -> Self: - """Create an AntibodyPanel from an AnnData object. + def concatenate(cls, panels: Sequence[PNAAntibodyPanel]) -> PNAAntibodyPanel: + """Concatenate panels into one panel. + + One panel is returned unchanged. Several panels are stacked in the + given order. The result keeps every input source. + ``marker_id``, ``sequence_1``, and ``sequence_2`` must be unique + across the concatenation. An optional column present on only some + sources is blank on the others, the same as an empty cell in a panel CSV. + """ + if not panels: + raise ValueError("At least one panel is required to concatenate.") + if len(panels) == 1: + return panels[0] + + frames: list[pd.DataFrame] = [] + sources: list[PanelSource] = [] + source_id_frames: list[pd.Series] = [] + for panel in panels: + if not panel.sources: + raise ValueError("Cannot concatenate a panel that has no sources.") + for source_index, source in enumerate(panel.sources): + new_source_id = len(sources) + sources.append(source) + marker_index = panel.marker_source_ids.index[ + panel.marker_source_ids == source_index + ] + part = panel.df.loc[list(marker_index)] + frames.append(part) + source_id_frames.append( + pd.Series(new_source_id, index=part.index, dtype="int64") + ) - Args: - adata: An AnnData object containing panel information. - file_name: The optional name of the file from which the AnnData object was loaded. + df = pd.concat(cls._align_optional_columns(frames)) + df.index.name = cls._INDEX_COLUMN + if "control" in df.columns: + df["control"] = df["control"].map( + lambda value: bool(value) if pd.notna(value) else False + ) + marker_source_ids = pd.concat(source_id_frames) + marker_source_ids.index = df.index + return cls( + df, + sources=sources, + marker_source_ids=marker_source_ids.astype("int64"), + ) - Returns: - The AntibodyPanel object. (AntibodyPanel) + @staticmethod + def _align_optional_columns(frames: list[pd.DataFrame]) -> list[pd.DataFrame]: + """Give every frame the same columns before they are stacked. - Raises: - KeyError: exception if panel information is missing in the AnnData object. + ``pd.concat`` inserts NaN where a column exists on only some frames. + A blank UniProt id is an empty string, and a missing hashing flag is + false. Filling after the stack mixes those types and breaks validation. """ - logger.debug("Creating Antibody panel from AnnData object") - try: - panel_metadata = adata.uns["panel_metadata"] - except KeyError as err: - logger.error( # pylint: disable=logging-not-lazy - f"The provided AnnData object does not contain {err}. " - + "Please, regenerate your data with the most recent version of pixelator." - ) - raise - panel_columns = panel_metadata.get("panel_columns") - if not panel_columns: - raise KeyError( - "The provided AnnData object does not contain panel columns information in the metadata. " - + "Please, regenerate your data with the most recent version of pixelator." - ) - df = adata.var[panel_columns] - metadata = AntibodyPanelMetadata.model_validate(panel_metadata) + columns = list( + dict.fromkeys(column for frame in frames for column in frame.columns) + ) + missing_on_some = [ + column + for column in columns + if any(column not in frame.columns for frame in frames) + ] + if not missing_on_some: + return frames + fills = { + column: PNAAntibodyPanel._missing_optional_value(frames, column) + for column in missing_on_some + } + aligned: list[pd.DataFrame] = [] + for frame in frames: + missing = { + column: fill + for column, fill in fills.items() + if column not in frame.columns + } + if not missing: + aligned.append(frame) + continue + extra = pd.DataFrame(missing, index=frame.index) + aligned.append(pd.concat([frame, extra], axis=1)) + return aligned - logger.debug("Antibody panel from AnnData object created") - return cls(df, metadata, file_name=file_name) + @staticmethod + def _missing_optional_value(frames: list[pd.DataFrame], column: str): + """Return the blank value for a column some sources do not have.""" + if column == "control": + return False + present = [ + frame[column].dropna() for frame in frames if column in frame.columns + ] + values = pd.concat(present) if present else pd.Series(dtype=object) + if values.empty: + return "" + if pd.api.types.is_bool_dtype(values) or all( + isinstance(value, bool) for value in values + ): + return False + if pd.api.types.is_numeric_dtype(values): + return 0 + return "" + + def _single_metadata(self) -> AntibodyPanelMetadata | None: + """Return metadata when this panel has exactly one source.""" + if len(self.sources) != 1: + return None + return self.sources[0].metadata + + def _require_single_metadata(self, field: str) -> AntibodyPanelMetadata: + """Return the only source metadata, refusing a concatenated panel.""" + metadata = self._single_metadata() + if metadata is None: + raise ValueError( + f"Panel {field} is only available for a single source. " + "Read it from each entry in sources." + ) + return metadata + + @property + def metadata(self) -> AntibodyPanelMetadata: + """Metadata for the only source in this panel. + + Raises: + ValueError: When this panel does not have exactly one source. + Read ``sources`` instead. + """ + return self._require_single_metadata("metadata") @property def name(self) -> str: - """Panel name from metadata. + """Panel name. - Returns: - The panel name. + Raises: + ValueError: When this panel does not have exactly one source. + Read ``sources`` instead. """ - return self.metadata.name + return self._require_single_metadata("name").name @property def product(self) -> Optional[str]: """Product identifier from metadata, if present. Returns: - Product name, or None when not provided in panel metadata. + Product name, or None when the single source does not set one. + + Raises: + ValueError: When this panel does not have exactly one source. + Read ``sources`` instead. """ - return self.metadata.product + return self._require_single_metadata("product").product @property def version(self) -> str: - """Panel version from metadata. + """Panel version. - Returns: - Semantic version string for this panel. + Raises: + ValueError: When this panel does not have exactly one source. + Read ``sources`` instead. """ - return self.metadata.version + return self._require_single_metadata("version").version @property def description(self) -> Optional[str]: - """Return the panel file description.""" - return self.metadata.description + """Return the panel file description. + + Raises: + ValueError: When this panel does not have exactly one source. + Read ``sources`` instead. + """ + return self._require_single_metadata("description").description @property def aliases(self) -> list[str]: - """Return the (optional) list of panel file aliases.""" - return self.metadata.aliases + """Return the (optional) list of panel file aliases. + + Raises: + ValueError: When this panel does not have exactly one source. + Read ``sources`` instead. + """ + return self._require_single_metadata("aliases").aliases @property def archived(self) -> Optional[bool]: - """Return whether the panel is marked as archived.""" - return self.metadata.archived + """Return whether the panel is marked as archived. + + Raises: + ValueError: When this panel does not have exactly one source. + Read ``sources`` instead. + """ + return self._require_single_metadata("archived").archived + + @property + def marker_source_ids(self) -> pd.Series: + """Return the panel-source index of each marker.""" + return self._marker_source_ids + + @property + def hashing_marker_ids(self) -> set[str]: + """Return marker ids flagged by the ``sample_hashing`` column.""" + return _hashing_marker_ids(self.df) @classmethod def _parse_header(cls, file: Path) -> AntibodyPanelMetadata: @@ -227,6 +466,7 @@ def _parse_header(cls, file: Path) -> AntibodyPanelMetadata: @classmethod def _parse_panel(cls, panel_file: Path) -> pd.DataFrame: + """Read the marker table from a panel CSV and convert ``control`` to bool.""" panel = pd.read_csv(str(panel_file), comment="#", index_col="marker_id").fillna( "" ) @@ -260,6 +500,63 @@ def filepath(self) -> Optional[Path]: """Return the full path of the marker panel file, if any.""" return self._filepath + def copy(self) -> Self: + """Return a shallow copy with its own marker table and source list.""" + return type(self)( + self.df.copy(), + file_name=self.filename, + filepath=self.filepath, + sources=list(self.sources), + marker_source_ids=self.marker_source_ids.copy(), + ) + + def source_as_panel(self, source_index: int) -> Self: + """Return one source as its own single-source panel.""" + source = self.sources[source_index] + marker_index = self.marker_source_ids.index[ + self.marker_source_ids == source_index + ] + df = self.df.loc[list(marker_index)].copy() + return type(self)( + df, + source.metadata, + file_name=source.file_name, + filepath=source.filepath, + ) + + def replace_source(self, source_index: int, replacement: PNAAntibodyPanel) -> Self: + """Return a copy with one source replaced by a single-source panel. + + An optional column present on only one side is blank on the other, + the same as an empty cell in a panel CSV. + """ + if len(replacement.sources) != 1: + raise ValueError("Replacement panel must come from a single source.") + keep = self.marker_source_ids.index[self.marker_source_ids != source_index] + kept = self.df.loc[list(keep)] + if kept.empty: + df = replacement.df.copy() + else: + kept, incoming = self._align_optional_columns([kept, replacement.df]) + df = pd.concat([kept, incoming]) + df.index.name = self._INDEX_COLUMN + source_ids = pd.concat( + [ + self.marker_source_ids.loc[list(keep)], + pd.Series(source_index, index=replacement.df.index, dtype="int64"), + ] + ) + sources = list(self.sources) + sources[source_index] = replacement.sources[0] + single_source = len(sources) == 1 + return type(self)( + df, + file_name=self.filename if single_source else None, + filepath=self.filepath if single_source else None, + sources=sources, + marker_source_ids=source_ids.astype("int64"), + ) + @cached_property def size(self) -> int: """Return the size of the marker panel.""" @@ -267,6 +564,7 @@ def size(self) -> int: @staticmethod def _validate_sequences(panel_df, sequence_col): + """Return errors when sequences differ in length or are not ATCG.""" errors = [] sequences = panel_df[sequence_col] ref_length = len(sequences.iloc[0]) @@ -283,6 +581,7 @@ def _validate_sequences(panel_df, sequence_col): @staticmethod def _validate_marker_names(panel_df): + """Return errors when marker ids contain underscores or whitespace.""" errors = [] if any(panel_df.index.str.contains("_")): # Markers should not contain underscores since this messes @@ -349,6 +648,13 @@ def validate_antibody_panel( errors.append(f"`{cls._INDEX_COLUMN}` is missing or is not set as index") return errors + if panel_df.index.duplicated().any(): + duplicated = panel_df.index[panel_df.index.duplicated()].unique().tolist() + errors.append( + "All values in column: marker_id were not unique. " + f"Offending values: {duplicated}" + ) + errors += cls._validate_marker_names(panel_df) if panel_df["control"].dtype != bool: @@ -380,22 +686,81 @@ def check_id(id_str): errors += cls._validate_sequences(panel_df, "sequence_1") errors += cls._validate_sequences(panel_df, "sequence_2") + errors += cls._validate_hashing_marker_ids(panel_df) + + return errors + @staticmethod + def _validate_hashing_marker_ids(panel_df: pd.DataFrame) -> list[str]: + """Return errors when hashing ids lack a ``-`` suffix or nest.""" + hashing_ids = _hashing_marker_ids(panel_df) + if not hashing_ids: + return [] + errors: list[str] = [] + missing_suffix = [ + marker_id + for marker_id in sorted(hashing_ids) + if split_hashing_marker_id(marker_id) is None + ] + if missing_suffix: + errors.append( + "Hashing marker ids must end with - (e.g. B2M-1). " + f"Offending values: {missing_suffix}" + ) + nested = sorted( + marker_id + for marker_id in hashing_ids + if (base := collapsed_hashing_marker_id(marker_id)) != marker_id + and base in hashing_ids + ) + if nested: + errors.append( + "Hashing marker ids must not collapse to another hashing id " + f"(e.g. B2M-1-1 next to B2M-1). Offending values: {nested}" + ) return errors def to_polars(self) -> pl.DataFrame: """Convert the panel to a Polars DataFrame.""" return pl.from_pandas(self.df, include_index=True) + def _source_marker_groups( + self, + ) -> list[tuple[AntibodyPanelMetadata, frozenset[str]]]: + """Return each source with its marker ids, independent of source order.""" + groups = [] + for source_index, source in enumerate(self.sources): + marker_ids = frozenset( + str(marker_id) + for marker_id in self.marker_source_ids.index[ + self.marker_source_ids == source_index + ] + ) + groups.append((source.metadata, marker_ids)) + return sorted( + groups, + key=lambda group: (group[0].model_dump_json(), tuple(sorted(group[1]))), + ) + def __eq__(self, other: object) -> bool: - """Check if two panels are equal based on their dataframes and metadata. + """Return whether two panels describe the same sources and markers. + + Row order, column order, source order, and the file a source was + loaded from are ignored. Each marker must still belong to the same + source. Args: other: Panel to compare for equality. """ if not isinstance(other, PNAAntibodyPanel): raise ValueError("Can only compare with another PNAAntibodyPanel") - return self.df.equals(other.df) and self.metadata == other.metadata + if self._source_marker_groups() != other._source_marker_groups(): + return False + left = self.df.sort_index() + right = other.df.sort_index() + if set(left.columns) != set(right.columns): + return False + return left.equals(right[list(left.columns)]) def load_antibody_panel(config: PNAConfig, panel: PathType) -> PNAAntibodyPanel: @@ -418,6 +783,20 @@ def load_antibody_panel(config: PNAConfig, panel: PathType) -> PNAAntibodyPanel: return panel_obj +def load_antibody_panels( + config: PNAConfig, panels: PathType | Sequence[PathType] +) -> PNAAntibodyPanel: + """Load one or more panels and concatenate them. + + A single path or name is returned as that panel. Several inputs are + concatenated in order. + """ + if isinstance(panels, (str, os.PathLike)): + return load_antibody_panel(config, panels) + loaded = [load_antibody_panel(config, panel) for panel in panels] + return PNAAntibodyPanel.concatenate(loaded) + + class PNAAntibodyPanelDiff: """Class representing the differences between two PNAAntibodyPanel objects.""" @@ -429,7 +808,15 @@ def __init__(self, panel_1: PNAAntibodyPanel, panel_2: PNAAntibodyPanel) -> None Args: panel_1: The first panel to compare. panel_2: The second panel to compare. + + Raises: + ValueError: When either panel does not have exactly one source. """ + if len(panel_1.sources) != 1 or len(panel_2.sources) != 1: + raise ValueError( + "PNAAntibodyPanelDiff only compares panels with a single source. " + "Split a concatenated panel with source_as_panel first." + ) self.panel_1 = panel_1 self.panel_2 = panel_2 @@ -567,57 +954,311 @@ def removed_clones(self) -> pl.DataFrame: ] ) - def upgrade_adata(self, adata: AnnData) -> AnnData: - """Upgrade an AnnData object with the changes between the two panels. + def changed_marker_ids(self) -> dict[str, str]: + """Return ``old marker_id -> new marker_id`` for markers whose id changed.""" + if ( + "marker_id" not in self.joined.columns + or "marker_id_panel_2" not in self.joined.columns + ): + return {} + both = self.joined.filter( + pl.col("marker_id").is_not_null() + & pl.col("marker_id_panel_2").is_not_null() + ) + mapping: dict[str, str] = {} + for old, new in zip( + both["marker_id"].to_list(), both["marker_id_panel_2"].to_list() + ): + if str(old) != str(new): + mapping[str(old)] = str(new) + return mapping + + +def sample_calling_hashing_collapsed( + hashing_ids: set[str], + *, + adata: AnnData | None = None, + pxl_file_metadata: dict | None = None, +) -> bool: + """Return whether sample calling has already collapsed hashing clones. + + Explicit ``hashing_collapsed`` on the pixel file metadata wins. Older + files are inferred from ``original_hash_counts_*`` columns or from hashing + clones that are absent from ``var``. Without an AnnData and without that + key, there is nothing to infer from, so this returns False. + """ + if pxl_file_metadata is not None and "hashing_collapsed" in pxl_file_metadata: + return bool(pxl_file_metadata["hashing_collapsed"]) + if adata is None: + return False + if any(str(col).startswith("original_hash_counts_") for col in adata.obs.columns): + return True + if not hashing_ids: + return False + var_names = {str(name) for name in adata.var_names} + return hashing_ids.isdisjoint(var_names) + + +def align_panel_patches( + panels: list[PNAAntibodyPanel], + adatas: list[AnnData] | None = None, + *, + pxl_file_metadata: list[dict] | None = None, +) -> tuple[list[PNAAntibodyPanel], list[dict[str, str]], list[dict[str, str]]]: + """Bump each source to the newest patch carried by another panel. + + Sources match on name and product, and only when major and minor versions + agree. A source that is not already on a panel is left alone. + ``pxl_file_metadata`` is the pixel file metadata for each panel, in the same + order, and carries ``hashing_collapsed`` when sample calling wrote the file. - Args: - adata: An AnnData object containing panel information. - """ - adata_panel = PNAAntibodyPanel.from_adata(adata) - if self.panel_1 != adata_panel: - raise ValueError( - "The provided AnnData object does not match the panel. Cannot upgrade." - f"Expected panel {self.panel_2.name} v{self.panel_2.version}, but got panel {adata_panel.name} v{adata_panel.version}." + Returns: + Updated panel copies, marker renames to apply to stored data (var, + edgelist, proximity, layouts), and hashing-clone renames for + ``original_hash_counts_*`` columns. Renames are keyed by input order. + """ + updated = [panel.copy() for panel in panels] + data_renames: list[dict[str, str]] = [{} for _ in panels] + hash_renames: list[dict[str, str]] = [{} for _ in panels] + + families: dict[tuple[str, str, tuple[int, ...]], list[tuple[int, int, Version]]] = ( + defaultdict(list) + ) + for panel_index, panel in enumerate(updated): + for source_index, source in enumerate(panel.sources): + identity = _source_identity(source) + if identity is None: + continue + version = Version(source.metadata.version) + minor = version.release[:2] + families[(*identity, minor)].append((panel_index, source_index, version)) + + for members in families.values(): + latest_idx, latest_source_idx, latest_version = max( + members, key=lambda member: member[2] + ) + latest_panel = updated[latest_idx].source_as_panel(latest_source_idx) + for panel_index, source_index, version in members: + if version == latest_version: + continue + current = updated[panel_index].source_as_panel(source_index) + logger.info( + "Upgrading panel source %s %s from %s to %s.", + current.name, + current.product, + current.version, + latest_panel.version, ) - - non_panel_columns = adata.var.copy()[ - [ - col - for col in adata.var.columns - if col not in adata.uns["panel_metadata"]["panel_columns"] - ] - + self.join_on_columns - ] - adata.var = ( - self.joined.select( - list( - set( - self.join_on_columns - + self.identical_columns - + [f"{col}_panel_2" for col in self.changed_columns] - + self.added_columns - ) + diff = PNAAntibodyPanelDiff(current, latest_panel) + clone_map = diff.changed_marker_ids() + _validate_hashing_renames(current, latest_panel, clone_map) + hashing_ids = current.hashing_marker_ids + pxl_metadata = ( + None if pxl_file_metadata is None else pxl_file_metadata[panel_index] + ) + collapsed = sample_calling_hashing_collapsed( + hashing_ids, + adata=None if adatas is None else adatas[panel_index], + pxl_file_metadata=pxl_metadata, + ) + if adatas is not None: + _require_expected_markers( + adata=adatas[panel_index], + old_panel=current, + new_panel=latest_panel, + clone_map=clone_map, + collapsed=collapsed, ) + data_map = _data_rename_map(clone_map, hashing_ids, collapsed=collapsed) + _merge_renames(data_renames[panel_index], data_map) + _merge_renames( + hash_renames[panel_index], + {old: new for old, new in clone_map.items() if old in hashing_ids}, ) - .rename({f"{col}_panel_2": col for col in self.changed_columns}) - # keep order and append new to the end - .select( - ["marker_id"] # index not in panel_metadata panel_columns below - + adata.uns["panel_metadata"]["panel_columns"] - + self.added_columns + updated[panel_index] = updated[panel_index].replace_source( + source_index, latest_panel ) - .to_pandas() - .set_index("marker_id") - ) - if adata.var.shape[0] != non_panel_columns.shape[0]: + + return updated, data_renames, hash_renames + + +def aligned_dataset_panel(panels: list[PNAAntibodyPanel]) -> PNAAntibodyPanel: + """Return one panel after per-source patch alignment across files. + + Several files that describe the same sources collapse to a single panel. + The order of ``--panel`` inputs does not have to match. + """ + if len(panels) == 1: + return panels[0] + updated, _, _ = align_panel_patches(panels) + first = updated[0] + for other in updated[1:]: + if first != other: raise ValueError( - "Row count mismatch in automatic patch panel patch version bump." + "Samples do not share the same panel sources after patch alignment." ) - adata.var = adata.var.join( - non_panel_columns.set_index(self.join_on_columns), - how="outer", - on=self.join_on_columns, + return first + + +def _source_identity(source: PanelSource) -> tuple[str, str] | None: + """Return ``(name, product)`` for a source, or None when product is unset.""" + product = source.metadata.product + if not product: + return None + return (source.metadata.name, product) + + +def _data_rename_map( + clone_map: dict[str, str], hashing_ids: set[str], *, collapsed: bool +) -> dict[str, str]: + """Return marker renames, collapsing hashing clones when ``collapsed`` is set.""" + if not collapsed: + return dict(clone_map) + data_map = {old: new for old, new in clone_map.items() if old not in hashing_ids} + families: dict[str, str] = {} + for old, new in clone_map.items(): + if old not in hashing_ids: + continue + old_base = collapsed_hashing_marker_id(old) + new_base = collapsed_hashing_marker_id(new) + if old_base in families and families[old_base] != new_base: + raise ValueError( + f"Hashing markers with collapsed name {old_base!r} do not share " + "a single new base name." + ) + families[old_base] = new_base + for old_base, new_base in families.items(): + if old_base != new_base: + data_map[old_base] = new_base + return data_map + + +def _validate_hashing_renames( + old_panel: PNAAntibodyPanel, + new_panel: PNAAntibodyPanel, + clone_map: dict[str, str], +) -> None: + """Reject a hashing rename that breaks the collapsed-name rules. + + Each of these fails: + + * ``B2M-1`` → ``C2M``: a hashing id must end with ``-``. + * ``B2M-1`` → ``C2M-2``: the numeric suffix stays. + * ``B2M-1`` → ``C2M-1`` and ``B2M-2`` → ``D2M-2``: clones in one family + share one new base. + * ``B2M`` stays while ``B2M-1`` → ``C2M-1``, or ``B2M`` → ``C2M`` while + ``B2M-1`` stays: a non-hashing marker that already has the collapsed + name renames with the family. + * ``CD19`` → ``C2M`` while ``B2M-1`` → ``C2M-1``: the family must not + land on a different non-hashing marker. + """ + old_hashing = old_panel.hashing_marker_ids + new_hashing = new_panel.hashing_marker_ids + if not old_hashing: + return + + families: dict[str, str] = {} + for old in sorted(old_hashing): + new = clone_map.get(old, old) + if new not in new_hashing and old not in clone_map: + continue + old_parts = split_hashing_marker_id(old) + new_parts = split_hashing_marker_id(new) + if old_parts is None or new_parts is None: + raise ValueError( + "Hashing marker ids must end with - to be renamed " + f"({old!r} -> {new!r})." + ) + if old_parts[1] != new_parts[1]: + raise ValueError( + "Hashing marker rename may only change the base name, not the " + f"hash group suffix {old_parts[1]}. Got {old!r} -> {new!r}." + ) + old_base = collapsed_hashing_marker_id(old) + new_base = collapsed_hashing_marker_id(new) + if old_base in families and families[old_base] != new_base: + raise ValueError( + f"Hashing markers with collapsed name {old_base!r} must keep a " + "single base name per hash group family." + ) + families[old_base] = new_base + + non_hashing = { + str(marker_id) + for marker_id in old_panel.markers + if str(marker_id) not in old_hashing + } + for old_base, new_base in families.items(): + if old_base in non_hashing: + moved = clone_map.get(old_base, old_base) + if moved != new_base: + # TODO: we should think about if we want this check or not in the long run. + # adding it here now we can always remove it later if we decide it's too strict. + raise ValueError( + f"Hashing collapsed base {old_base!r} and non-hashing marker " + f"{old_base!r} must be renamed together " + f"(hashing maps to {new_base!r}, non-hashing maps to {moved!r})." + ) + for other in non_hashing: + if other == old_base: + continue + if clone_map.get(other, other) == new_base: + raise ValueError( + f"Hashing collapsed base {old_base!r} -> {new_base!r} collides " + f"with non-hashing marker {other!r}." + ) + + +def _require_expected_markers( + *, + adata: AnnData, + old_panel: PNAAntibodyPanel, + new_panel: PNAAntibodyPanel, + clone_map: dict[str, str], + collapsed: bool, +) -> None: + """Raise when a patch bump is missing a marker the file should still contain. + + After sample calling, hashing clones may be absent when the file is + collapsed. A missing non-hashing marker still fails. + """ + var_names = {str(name) for name in adata.var_names} + old_hashing = old_panel.hashing_marker_ids + new_hashing = new_panel.hashing_marker_ids + old_markers = {str(marker_id) for marker_id in old_panel.markers} + missing: list[str] = [] + for marker_id in old_panel.markers: + marker = str(marker_id) + if marker in var_names or marker in clone_map: + continue + if collapsed and marker in old_hashing: + continue + missing.append(marker) + for marker_id in new_panel.markers: + marker = str(marker_id) + if marker in var_names or marker in old_markers: + continue + old_ids = [old for old, new in clone_map.items() if new == marker] + if any(old in var_names for old in old_ids): + continue + if collapsed and ( + marker in new_hashing or any(old in old_hashing for old in old_ids) + ): + continue + missing.append(marker) + if missing: + raise ValueError( + "Row count mismatch in automatic panel patch version bump. " + f"Missing markers: {sorted(set(missing))[:5]}" ) - adata.uns["panel_metadata"] = self.panel_2.metadata.model_dump() - adata.uns["panel_metadata"]["panel_columns"] = self.panel_2.df.columns.tolist() - return adata + + +def _merge_renames(target: dict[str, str], extra: dict[str, str]) -> None: + """Copy renames into ``target``, rejecting two destinations for one marker.""" + for old, new in extra.items(): + if old in target and target[old] != new: + raise ValueError( + f"Marker {old!r} is renamed to both {target[old]!r} and {new!r}." + ) + if old != new: + target[old] = new diff --git a/src/pixelator/pna/config/panel_align.py b/src/pixelator/pna/config/panel_align.py new file mode 100644 index 000000000..6f97121da --- /dev/null +++ b/src/pixelator/pna/config/panel_align.py @@ -0,0 +1,45 @@ +"""Apply a panel patch bump to AnnData marker ids. + +Copyright © 2026 Pixelgen Technologies AB. +""" + +from __future__ import annotations + +from typing import TYPE_CHECKING + +if TYPE_CHECKING: + from anndata import AnnData + + +def apply_marker_renames_to_adata(adata: AnnData, mapping: dict[str, str]) -> None: + """Rename ``var`` marker ids in place.""" + if not mapping: + return + new_names = [mapping.get(str(name), str(name)) for name in adata.var_names] + if len(new_names) != len(set(new_names)): + raise ValueError( + "Panel patch bump would duplicate marker ids in var: " + f"{sorted(name for name in new_names if new_names.count(name) > 1)}" + ) + adata.var_names = new_names + + +def apply_hash_count_renames(adata: AnnData, mapping: dict[str, str]) -> None: + """Rename ``original_hash_counts_*`` columns when hashing ids change.""" + if not mapping: + return + renames = {} + for old, new in mapping.items(): + if old == new: + continue + old_col = f"original_hash_counts_{old}" + new_col = f"original_hash_counts_{new}" + if old_col in adata.obs.columns: + if new_col in adata.obs.columns: + raise ValueError( + f"Cannot rename hashing counts column {old_col} to {new_col}: " + "the target already exists." + ) + renames[old_col] = new_col + if renames: + adata.obs.rename(columns=renames, inplace=True) diff --git a/src/pixelator/pna/config/panel_tables.py b/src/pixelator/pna/config/panel_tables.py new file mode 100644 index 000000000..7cbe9afd7 --- /dev/null +++ b/src/pixelator/pna/config/panel_tables.py @@ -0,0 +1,97 @@ +"""DuckDB storage for panel markers and their sources. + +Copyright © 2026 Pixelgen Technologies AB. +""" + +from __future__ import annotations + +import json +from typing import TYPE_CHECKING + +import pandas as pd + +if TYPE_CHECKING: + import duckdb + + from pixelator.pna.config.panel import PNAAntibodyPanel + + +_PANEL_TABLES = ("panels", "panel_sources") +_PANEL_STORAGE_COLUMNS = frozenset({"row_nr", "source_id", "marker_id"}) + + +def panel_tables_present(connection: duckdb.DuckDBPyConnection) -> bool: + """Return whether both panel tables exist on this connection.""" + names = set(connection.execute("SHOW TABLES").fetchdf()["name"].tolist()) + return set(_PANEL_TABLES).issubset(names) + + +def stored_panel_marker_columns(connection: duckdb.DuckDBPyConnection) -> set[str]: + """Return marker columns stored in ``panels``. + + Empty when the panel tables are absent. ``row_nr``, ``source_id``, and + ``marker_id`` are storage columns, not marker fields. + """ + if not panel_tables_present(connection): + return set() + frame = connection.execute("SELECT * FROM panels LIMIT 0").fetchdf() + return {str(name) for name in frame.columns if name not in _PANEL_STORAGE_COLUMNS} + + +def write_panel_tables( + connection: duckdb.DuckDBPyConnection, panel: PNAAntibodyPanel +) -> None: + """Replace ``panels`` and ``panel_sources`` with ``panel``.""" + if not panel.sources: + raise ValueError("Cannot store a panel that has no sources.") + + source_rows = [] + for source_id, source in enumerate(panel.sources): + metadata = source.metadata + source_rows.append( + { + "source_id": source_id, + "name": metadata.name, + "version": metadata.version, + "product": metadata.product, + "description": metadata.description, + "aliases": json.dumps(list(metadata.aliases)), + "archived": bool(metadata.archived) + if metadata.archived is not None + else False, + "file_name": source.file_name, + "filepath": source.filepath, + } + ) + sources_df = pd.DataFrame(source_rows) + + markers = panel.df.copy() + markers.insert(0, "marker_id", markers.index.astype(str)) + markers.insert(0, "row_nr", range(len(markers))) + markers["source_id"] = ( + panel.marker_source_ids.reindex(panel.df.index).astype(int).to_numpy() + ) + + connection.register("_panel_sources_df", sources_df) + connection.register("_panels_df", markers.reset_index(drop=True)) + try: + connection.execute( + "CREATE OR REPLACE TABLE panel_sources AS SELECT * FROM _panel_sources_df" + ) + connection.execute("CREATE OR REPLACE TABLE panels AS SELECT * FROM _panels_df") + finally: + connection.unregister("_panel_sources_df") + connection.unregister("_panels_df") + + +def read_panel_table_frames( + connection: duckdb.DuckDBPyConnection, +) -> tuple[pd.DataFrame, pd.DataFrame] | None: + """Return ``(markers, sources)`` or None when the tables are absent.""" + if not panel_tables_present(connection): + return None + markers = connection.execute("SELECT * FROM panels ORDER BY row_nr").fetchdf() + sources = connection.execute( + "SELECT * FROM panel_sources ORDER BY source_id" + ).fetchdf() + return markers, sources diff --git a/src/pixelator/pna/graph/component_recovery.py b/src/pixelator/pna/graph/component_recovery.py index d4e756ddc..5a4137fae 100644 --- a/src/pixelator/pna/graph/component_recovery.py +++ b/src/pixelator/pna/graph/component_recovery.py @@ -122,13 +122,12 @@ def build_pxl_file_with_components( "sample_name": sample_name, "version": __version__, "technology": "single-cell-pna", - "panel_name": panel.name, - "panel_version": panel.version, } logger.debug("Building pxl file") pxl_file_writer.write_metadata(metadata) + pxl_file_writer.write_panel(panel) pxl_file_writer.write_adata(adata) return PNAPixelDataset.from_pxl_files(path_output_pxl_file), component_stats diff --git a/src/pixelator/pna/pixeldataset/edgelist.py b/src/pixelator/pna/pixeldataset/edgelist.py index 3bb7e0dab..ced4a8424 100644 --- a/src/pixelator/pna/pixeldataset/edgelist.py +++ b/src/pixelator/pna/pixeldataset/edgelist.py @@ -99,12 +99,33 @@ def to_polars(self) -> pl.DataFrame: def to_record_batches( self, batch_size: int = 1_000_000 ) -> Iterable[pa.RecordBatch]: - """Get the edgelist as a stream of pyarrow RecordBatches.""" + """Get the edgelist as a stream of pyarrow RecordBatches. + + Legacy ``marker1`` and ``marker2`` columns are renamed to + ``marker_1`` and ``marker_2``, matching :meth:`to_polars`. + """ query = self._query_builder.edgelist_query( normalize_input_to_list(self.components) ) with self._view.open() as session: - yield from session.execute_arrow_reader(query=query, batch_size=batch_size) + for batch in session.execute_arrow_reader( + query=query, batch_size=batch_size + ): + names = set(batch.schema.names) + if "marker1" not in names and "marker2" not in names: + yield batch + continue + frame = pl.from_arrow(batch) + if not isinstance(frame, pl.DataFrame): + yield batch + continue + renamed = self._handle_backwards_compatibility(frame.lazy()).collect() + table = renamed.to_arrow() + batches = table.to_batches(max_chunksize=max(batch.num_rows, 1)) + if batches: + yield from batches + else: + yield pa.RecordBatch.from_pylist([], schema=table.schema) def _iterator(self) -> Iterable[tuple[str, pl.LazyFrame]]: with self._view.open() as session: diff --git a/src/pixelator/pna/pixeldataset/io/__init__.py b/src/pixelator/pna/pixeldataset/io/__init__.py index e611b891e..9ae2c5aa0 100644 --- a/src/pixelator/pna/pixeldataset/io/__init__.py +++ b/src/pixelator/pna/pixeldataset/io/__init__.py @@ -119,6 +119,7 @@ PXL_FILE_OTHER_TABLES, NullPxlFileError, PxlFile, + read_dataset_panel, reject_null_pxl, write_null_pxl, ) @@ -135,6 +136,7 @@ "PxlFile", "Query", "QueryBuilder", + "read_dataset_panel", "reject_null_pxl", "write_null_pxl", ] diff --git a/src/pixelator/pna/pixeldataset/io/anndata_helper.py b/src/pixelator/pna/pixeldataset/io/anndata_helper.py index ee42a40b1..63de51770 100644 --- a/src/pixelator/pna/pixeldataset/io/anndata_helper.py +++ b/src/pixelator/pna/pixeldataset/io/anndata_helper.py @@ -11,13 +11,10 @@ from typing import Literal import polars as pl -import semver from anndata import AnnData, ImplicitModificationWarning from anndata import concat as anndata_concat from pixelator.common.statistics import clr_transformation, log1p_transformation -from pixelator.common.utils import logger -from pixelator.pna.config.panel import PNAAntibodyPanel, PNAAntibodyPanelDiff from pixelator.pna.pixeldataset.utils import update_metrics_anndata from pixelator.pna.utils import normalize_input_to_list, normalize_input_to_set @@ -54,13 +51,14 @@ def _read_all_samples(self) -> AnnData: """ with self._view.open() as session: adatas: list[AnnData] = [] - for sample_name in self._view.sample_names(): + samples = list(self._view.sample_names()) + for sample_name in samples: adata = self._read_adata_from_sample( session=session, sample=sample_name ) adatas.append(adata) - self._try_bump_adata_panel_version(adatas) + self._view.apply_panel_patch_to_adatas(samples, adatas) if not adatas: return AnnData() @@ -76,49 +74,6 @@ def _read_all_samples(self) -> AnnData: update_metrics_anndata(concatenated, inplace=True) return concatenated - def _try_bump_adata_panel_version( - self, - adatas: list[AnnData], - ) -> list[AnnData]: - """Try to bump the panel version of the given AnnData. - - Only try to upgrade to the latest version available in the view, - if the panels differ in patch version and have the same product. - """ - if any("panel_metadata" not in adata.uns for adata in adatas): - logger.debug( - "Missing panel metadata in one or more samples, " - + "skipping automatic panel version upgrade." - ) - return adatas - panel_versions = [ - semver.Version.parse(adata.uns["panel_metadata"]["version"]) - for adata in adatas - ] - panel_products = [ - adata.uns["panel_metadata"].get("product") for adata in adatas - ] - latest_panel_version = max(panel_versions) - if ( - all(ver.major == latest_panel_version.major for ver in panel_versions) - and all(ver.minor == latest_panel_version.minor for ver in panel_versions) - and len(set(panel_products)) == 1 - and panel_products[0] is not None - ): - logger.info( - "Multiple panel patch versions for same product detected across samples. " - + "Attempting to upgrade to the latest." - ) - latest_panel = PNAAntibodyPanel.from_adata( - adatas[panel_versions.index(latest_panel_version)] - ) - for i, adata in enumerate(adatas): - panel = PNAAntibodyPanel.from_adata(adata) - if panel.version != latest_panel.version: - diff = PNAAntibodyPanelDiff(panel, latest_panel) - adatas[i] = diff.upgrade_adata(adata) - return adatas - def _read_adata_from_sample( self, *, diff --git a/src/pixelator/pna/pixeldataset/io/pixel_data_viewer.py b/src/pixelator/pna/pixeldataset/io/pixel_data_viewer.py index 4368aec67..eaf9397b8 100644 --- a/src/pixelator/pna/pixeldataset/io/pixel_data_viewer.py +++ b/src/pixelator/pna/pixeldataset/io/pixel_data_viewer.py @@ -11,8 +11,14 @@ import duckdb import polars as pl +from anndata import AnnData from pixelator.common.duckdb_utils import connect_duckdb +from pixelator.pna.config.panel import PNAAntibodyPanel, align_panel_patches +from pixelator.pna.config.panel_align import ( + apply_hash_count_renames, + apply_marker_renames_to_adata, +) from .pxl_file import PXL_FILE_MANDATOR_TABLES, PXL_FILE_OTHER_TABLES, PxlFile from .query_builder import Query @@ -23,6 +29,31 @@ # Characters / substrings that must not appear in sample labels embedded as '...' literals. _SAMPLE_LABEL_UNSAFE_RE = re.compile(r"['\";\x00-\x1f\\]|--|/\*|\*/") +# Marker-id columns in the union views. A patch bump renames these in the +# view, so later queries see the bumped ids. The attached files are unchanged. +_TABLE_MARKER_COLUMNS: dict[str, tuple[str, ...]] = { + "edgelist": ("marker_1", "marker_2", "marker1", "marker2"), + "proximity": ("marker_1", "marker_2"), + "layouts": ("marker",), +} + + +def _sql_string_literal(value: str) -> str: + """Return a single-quoted SQL literal.""" + return "'" + value.replace("'", "''") + "'" + + +def _marker_rename_expr(column: str, old_to_new: dict[str, str]) -> str | None: + """Return a CASE expression that maps stored marker ids to bumped ids.""" + branches = [ + f"WHEN {_sql_string_literal(old)} THEN {_sql_string_literal(new)}" + for old, new in old_to_new.items() + if old != new + ] + if not branches: + return None + return f"CASE {column} {' '.join(branches)} ELSE {column} END" + def _validate_session_sql_identifier(value: str, *, what: str) -> None: if not _SESSION_SQL_IDENTIFIER_RE.fullmatch(value): @@ -89,6 +120,10 @@ def __init__( if invalid_files: raise ValueError(f"{invalid_files} are not valid PXL files.") + self._marker_renames: dict[str, dict[str, str]] | None = None + self._hash_count_renames: dict[str, dict[str, str]] | None = None + self._upgraded_panels: dict[str, PNAAntibodyPanel | None] | None = None + def _map_sample_names_to_db_names( self, sample_name_to_pxl_file_mapping: dict[str, PxlFile] ) -> dict[str, str]: @@ -148,13 +183,34 @@ def from_sample_to_file_mappings( return PixelDataViewer(sample_name_to_pxl_file_mapping) def filter_samples(self, sample_names: set[str]) -> "PixelDataViewer": - """Create a new PixelDataViewer with only a subset of the samples.""" + """Create a new PixelDataViewer with only a subset of the samples. + + When this view holds several files, the new view keeps the patch bump + settled for those files. A later read of one remaining sample still + uses those marker names. + """ filtered_mapping = { sample_name: pxl_file for sample_name, pxl_file in self._db_to_file_mapping.items() if sample_name in sample_names } - return PixelDataViewer.from_sample_to_file_mappings(filtered_mapping) + viewer = PixelDataViewer.from_sample_to_file_mappings(filtered_mapping) + if len(self._db_to_file_mapping) < 2: + return viewer + self.marker_renames_by_sample() + marker_renames = self._marker_renames or {} + hash_renames = self._hash_count_renames or {} + upgraded = self._upgraded_panels or {} + viewer._marker_renames = { + sample: marker_renames[sample] for sample in filtered_mapping + } + viewer._hash_count_renames = { + sample: hash_renames[sample] for sample in filtered_mapping + } + viewer._upgraded_panels = { + sample: upgraded[sample] for sample in filtered_mapping + } + return viewer @property def sample_to_file_mappings(self) -> dict[str, Path]: @@ -165,12 +221,27 @@ def sample_to_file_mappings(self) -> dict[str, Path]: } def open(self) -> PixelDataViewerSession: - """Return a new session with an open DuckDB connection (context manager or ``close()``).""" + """Return a session whose marker columns use the patch bump. + + Edgelist, proximity, and layouts are renamed in the session views. + The files on disk are not rewritten. A view opened on one file shows + the stored ids. A view narrowed from several files keeps that bump. + Computing the bump reads AnnData through a session that still shows + the stored ids. + """ + if self._marker_renames is None and len(self._db_to_file_mapping) >= 2: + self.marker_renames_by_sample() + return self._open_session(self._marker_renames or {}) + + def _open_session( + self, marker_renames: dict[str, dict[str, str]] | None = None + ) -> PixelDataViewerSession: + """Open a session, optionally renaming marker columns in its views.""" sources: list[tuple[str, Path, str]] = [ (sample_name, pxl_file.path, self._get_normalized_name(sample_name)) for sample_name, pxl_file in self._db_to_file_mapping.items() ] - return PixelDataViewerSession(sources) # type: ignore[arg-type] + return PixelDataViewerSession(sources, marker_renames=marker_renames) # type: ignore[arg-type] def sample_names(self) -> list[str]: """Return the list of sample names known to the view.""" @@ -180,6 +251,129 @@ def normalized_sample_db_name(self, sample_name: str) -> str: """Return the attached DuckDB database name for a sample.""" return self._get_normalized_name(sample_name) + def marker_renames_by_sample(self) -> dict[str, dict[str, str]]: + """Return old-to-new marker ids for each sample in this view. + + The map is the in-memory panel patch bump. Files on disk are unchanged. + """ + if self._marker_renames is None: + samples = list(self.sample_names()) + self._compute_panel_patch( + samples, self._read_adatas_for_panel_patch(samples) + ) + return dict(self._marker_renames or {}) + + def apply_panel_patch_to_adatas( + self, + samples: list[str] | list[AnnData], + adatas: list[AnnData] | None = None, + ) -> list[AnnData]: + """Apply the view's panel patch bump and attach panel columns. + + Marker ids are renamed first. Panel columns already on ``var`` are + then replaced by the panel this view settled on, so a legacy file and + a new file end up with the same in-memory columns. ``uns`` panel + metadata is removed. ``samples`` may be omitted by passing only the + AnnData list. Sample names are then taken from the view, in view order. + """ + if adatas is None: + adatas = list(samples) # type: ignore[arg-type] + samples = list(self.sample_names()) + sample_names = [str(sample) for sample in samples] + if self._marker_renames is None: + self._compute_panel_patch(sample_names, adatas) + else: + self._apply_cached_panel_patch(sample_names, adatas) + return adatas + + def _read_adatas_for_panel_patch(self, samples: list[str]) -> list[AnnData]: + # Imported here because AnnDataHelper imports this module. + from pixelator.pna.pixeldataset.io.anndata_helper import AnnDataHelper + + helper = AnnDataHelper(self) + with self._open_session() as session: + return [ + helper._read_adata_from_sample(session=session, sample=sample) + for sample in samples + ] + + def _compute_panel_patch(self, samples: list[str], adatas: list[AnnData]) -> None: + panels: list[PNAAntibodyPanel] = [] + indexed_adatas: list[AnnData] = [] + indexed_samples: list[str] = [] + indexed_metadata: list[dict] = [] + for sample, adata in zip(samples, adatas): + pxl_file = self._db_to_file_mapping[sample] + panel = pxl_file.read_panel() + if panel is None: + continue + panels.append(panel) + indexed_adatas.append(adata) + indexed_samples.append(sample) + indexed_metadata.append(pxl_file.metadata()) + + marker_renames: dict[str, dict[str, str]] = {sample: {} for sample in samples} + hash_count_renames: dict[str, dict[str, str]] = { + sample: {} for sample in samples + } + upgraded_panels: dict[str, PNAAntibodyPanel | None] = { + sample: None for sample in samples + } + for sample, panel in zip(indexed_samples, panels): + upgraded_panels[sample] = panel + if len(panels) >= 2: + upgraded, data_renames, hash_renames = align_panel_patches( + panels, indexed_adatas, pxl_file_metadata=indexed_metadata + ) + for sample, panel, data_map, hash_map in zip( + indexed_samples, upgraded, data_renames, hash_renames + ): + marker_renames[sample] = data_map + hash_count_renames[sample] = hash_map + upgraded_panels[sample] = panel + self._marker_renames = marker_renames + self._hash_count_renames = hash_count_renames + self._upgraded_panels = upgraded_panels + self._apply_cached_panel_patch(samples, adatas) + + def _apply_cached_panel_patch( + self, samples: list[str], adatas: list[AnnData] + ) -> None: + marker_renames = self._marker_renames or {} + hash_renames = self._hash_count_renames or {} + upgraded = self._upgraded_panels or {} + for sample, adata in zip(samples, adatas): + apply_marker_renames_to_adata(adata, marker_renames.get(sample, {})) + apply_hash_count_renames(adata, hash_renames.get(sample, {})) + panel = upgraded.get(sample) + if panel is not None: + _attach_panel_columns(adata, panel) + + +def _attach_panel_columns(adata: AnnData, panel: PNAAntibodyPanel) -> None: + """Join ``panel`` onto ``var`` and drop legacy ``uns`` panel metadata. + + Columns already on ``var`` that also belong to the panel are removed + first, so a legacy file does not keep the pre-bump values. A column the + stored panel listed and this panel no longer has is removed as well. + Markers that are not already in ``var`` are not added. + """ + metadata = adata.uns.get("panel_metadata") if adata.uns is not None else None + legacy_columns = set() + if isinstance(metadata, dict): + legacy_columns = {str(column) for column in metadata.get("panel_columns", [])} + panel_columns = set(map(str, panel.df.columns)) + to_drop = [ + column + for column in adata.var.columns + if str(column) in panel_columns or str(column) in legacy_columns + ] + if to_drop: + adata.var = adata.var.drop(columns=to_drop) + adata.var = adata.var.join(panel.df.reindex(adata.var_names), how="left") + if adata.uns is not None and "panel_metadata" in adata.uns: + del adata.uns["panel_metadata"] + class PixelDataViewerSession: """DuckDB session over one or more attached PXL files. @@ -195,8 +389,17 @@ class PixelDataViewerSession: calls ``close()``). Each session uses its own connection. """ - def __init__(self, sources: list[tuple[str, Path | str, str]]) -> None: - """Open a DuckDB connection and attach each PXL file in ``sources``.""" + def __init__( + self, + sources: list[tuple[str, Path | str, str]], + marker_renames: dict[str, dict[str, str]] | None = None, + ) -> None: + """Open a DuckDB connection and attach each PXL file in ``sources``. + + ``marker_renames`` maps each sample to its old-to-new marker ids. + The union views apply that map. An empty map leaves the stored ids. + """ + self._marker_renames = marker_renames or {} normalized: list[tuple[str, Path, str]] = [ (sample_name, Path(path), db_name) for sample_name, path, db_name in sources ] @@ -223,6 +426,31 @@ def _attach_to_files(self, connection: duckdb.DuckDBPyConnection) -> None: query += f"ATTACH DATABASE '{path}' AS {db_name} (READ_ONLY);\n" connection.execute(query) + def _sample_select( + self, + sample_name: str, + db_name: str, + table_name: str, + columns: set[str], + ) -> str: + """Select one sample's table, renaming marker ids when this session has a map.""" + mapping = self._marker_renames.get(sample_name, {}) + replacements = [] + for column in _TABLE_MARKER_COLUMNS.get(table_name, ()): + if column not in columns: + continue + expr = _marker_rename_expr(column, mapping) + if expr is not None: + replacements.append(f"{expr} AS {column}") + sample_sql = _sql_string_literal(sample_name) + if not replacements: + return f"SELECT *, {sample_sql} AS sample FROM {db_name}.{table_name}" + replaced = ", ".join(replacements) + return ( + f"SELECT * REPLACE ({replaced}), {sample_sql} AS sample " + f"FROM {db_name}.{table_name}" + ) + def _simple_union_table_view( self, connection: duckdb.DuckDBPyConnection, @@ -238,8 +466,12 @@ def _simple_union_table_view( # see: https://github.com/duckdb/duckdb/issues/13069 table_queries: list[str] = [] for sample_name, _path, db_name in self._sources: + described = connection.execute( + f"DESCRIBE {db_name}.{table_name}" + ).fetchall() + columns = {str(row[0]) for row in described} table_queries.append( - f"SELECT *, '{sample_name}' AS sample FROM {db_name}.{table_name}" + self._sample_select(sample_name, db_name, table_name, columns) ) if not table_queries: return diff --git a/src/pixelator/pna/pixeldataset/io/pixel_file_writer.py b/src/pixelator/pna/pixeldataset/io/pixel_file_writer.py index 48356ae14..9c671456f 100644 --- a/src/pixelator/pna/pixeldataset/io/pixel_file_writer.py +++ b/src/pixelator/pna/pixeldataset/io/pixel_file_writer.py @@ -5,6 +5,7 @@ from __future__ import annotations +import json import tempfile from pathlib import Path @@ -12,6 +13,12 @@ import polars as pl from anndata import AnnData +from pixelator.pna.config.panel import PNAAntibodyPanel +from pixelator.pna.config.panel_tables import ( + panel_tables_present, + stored_panel_marker_columns, + write_panel_tables, +) from pixelator.pna.utils import init_duckdb_conn @@ -173,6 +180,12 @@ def _clean_existing_adata_tables(self): def write_adata(self, adata: AnnData) -> None: """Write the AnnData object to the PXL file. + When the panel tables are already on this file, marker columns from + that panel are left out of ``var``. Those columns are joined back + on read. A file from pixelator 0.22.0 through 0.30.0 stores the + panel in ``uns`` instead. That panel is copied into the panel + tables before the AnnData tables are replaced. + Args: adata: The AnnData object to write. @@ -187,10 +200,16 @@ def write_adata(self, adata: AnnData) -> None: "but a PXL file may only contain a single sample." ) + self._store_legacy_panel() self._clean_existing_adata_tables() X = adata.to_df().reset_index(names="index") var = adata.var.reset_index(names="index") + panel_columns = stored_panel_marker_columns(self._connection) + if panel_columns: + var = var.drop( + columns=[column for column in var.columns if column in panel_columns] + ) obs = adata.obs.reset_index(names="index") uns = adata.uns @@ -214,6 +233,44 @@ def write_adata(self, adata: AnnData) -> None: """, ) + def _store_legacy_panel(self) -> None: + """Copy a 0.22.0–0.30.0 panel into the panel tables. + + The panel still sits in ``__adata__uns`` until this write replaces + those tables. Files that already have panel tables, and files from + before 0.22.0, are left unchanged. + """ + if panel_tables_present(self._connection): + return + try: + uns_row = self._connection.execute( + "SELECT value FROM __adata__uns" + ).fetchone() + except duckdb.CatalogException: + return + if uns_row is None: + return + uns = json.loads(uns_row[0]) if isinstance(uns_row[0], str) else uns_row[0] + if not isinstance(uns, dict) or "panel_metadata" not in uns: + return + try: + var = self._connection.execute("SELECT * FROM __adata__var").fetchdf() + except duckdb.CatalogException: + return + var = var.set_index("index").rename_axis(index={"index": "marker_id"}) + panel = PNAAntibodyPanel.from_legacy_var( + var, uns["panel_metadata"], file_name=self.path.name + ) + write_panel_tables(self._connection, panel) + + def write_panel(self, panel) -> None: + """Write ``panels`` and ``panel_sources`` for ``panel``. + + Args: + panel: The antibody panel to store. + """ + write_panel_tables(self.get_connection(), panel) + def write_metadata(self, metadata: dict) -> None: """Write the metadata to the PXL file. diff --git a/src/pixelator/pna/pixeldataset/io/pxl_file.py b/src/pixelator/pna/pixeldataset/io/pxl_file.py index 42d3c5bd8..70979dc7c 100644 --- a/src/pixelator/pna/pixeldataset/io/pxl_file.py +++ b/src/pixelator/pna/pixeldataset/io/pxl_file.py @@ -13,11 +13,17 @@ import json import shutil from pathlib import Path +from typing import TYPE_CHECKING import duckdb from pixelator.common.duckdb_utils import connect_duckdb from pixelator.common.exceptions import PixelatorBaseException +from pixelator.pna.config.panel import PNAAntibodyPanel, aligned_dataset_panel +from pixelator.pna.config.panel_tables import read_panel_table_frames + +if TYPE_CHECKING: + from pixelator.pna.pixeldataset.dataset import PNAPixelDataset PXL_FILE_MANDATOR_TABLES = [ "__adata__X", @@ -99,6 +105,35 @@ def metadata(self) -> dict: except duckdb.CatalogException: return {} + def read_panel(self) -> PNAAntibodyPanel | None: + """Load the panel stored in this file. + + New files are read from the ``panels`` and ``panel_sources`` tables. + Files from pixelator 0.22.0 through 0.30.0 are read from + ``uns['panel_metadata']`` and the panel columns stored on ``var``. + Returns None when the file has neither. + """ + with connect_duckdb(self.path, read_only=True) as connection: + frames = read_panel_table_frames(connection) + if frames is not None: + return PNAAntibodyPanel._from_panel_tables(*frames) + try: + uns_row = connection.execute( + "SELECT value FROM __adata__uns" + ).fetchone() + except duckdb.CatalogException: + return None + if uns_row is None: + return None + uns = json.loads(uns_row[0]) if isinstance(uns_row[0], str) else uns_row[0] + if not isinstance(uns, dict) or "panel_metadata" not in uns: + return None + var = connection.execute("SELECT * FROM __adata__var").fetchdf() + var = var.set_index("index").rename_axis(index={"index": "marker_id"}) + return PNAAntibodyPanel.from_legacy_var( + var, uns["panel_metadata"], file_name=self.path.name + ) + def is_null_file(self) -> bool: """Return True when this file is a null pxl file.""" return self.metadata().get("null") is True @@ -140,6 +175,32 @@ def copy_pxl_file(src: PxlFile, target: Path) -> PxlFile: return PxlFile(target) +def read_dataset_panel(dataset: PNAPixelDataset) -> PNAAntibodyPanel: + """Load one panel from every file in ``dataset``. + + One file is returned as that panel. Several files are aligned to the + newest patch of each source they share. + + Raises: + ValueError: If a file has no panel. Pixelator 0.22.0 through 0.30.0 + stored it in ``uns['panel_metadata']``. Later files store it in + the ``panels`` and ``panel_sources`` tables. Earlier files have + neither and must be rerun with a current version of the software. + """ + panels = [] + for path in dataset.view.sample_to_file_mappings.values(): + panel = PxlFile(path).read_panel() + if panel is None: + raise ValueError( + f"{path} was written before pixelator 0.22.0 and has no panel. " + "Rerun this sample with a current version of the software." + ) + panels.append(panel) + if len(panels) == 1: + return panels[0] + return aligned_dataset_panel(panels) + + _EMPTY_NULL_TABLES_SQL = """ CREATE TABLE edgelist ( umi1 UBIGINT, diff --git a/src/pixelator/pna/sample_calling/__init__.py b/src/pixelator/pna/sample_calling/__init__.py index fdcffc229..a8a4019a8 100644 --- a/src/pixelator/pna/sample_calling/__init__.py +++ b/src/pixelator/pna/sample_calling/__init__.py @@ -4,7 +4,6 @@ """ import logging -import re import tempfile from itertools import chain from pathlib import Path @@ -20,9 +19,9 @@ from pixelator import __version__ from pixelator.pna.analysis_engine import AnalysisManager, PerComponentTask from pixelator.pna.anndata import add_missing_adata_info, pna_edgelist_to_anndata -from pixelator.pna.config.panel import PNAAntibodyPanel +from pixelator.pna.config.panel import PNAAntibodyPanel, collapsed_hashing_marker_id from pixelator.pna.pixeldataset import PNAPixelDataset -from pixelator.pna.pixeldataset.io import PixelFileWriter +from pixelator.pna.pixeldataset.io import PixelFileWriter, read_dataset_panel from pixelator.pna.sample_calling.hash_antibodies import HashedAntibodyMapping from pixelator.pna.sample_calling.report import SampleCallingTotalReport @@ -209,6 +208,25 @@ def _add_original_hash_counts_to_obs( old_adata.obs[f"original_hash_counts_{ab}"] = 0 +def _count_markers_after_hash_collapse( + panel: PNAAntibodyPanel, hashing_antibodies: set[str] +) -> list[str]: + """Return count-matrix markers after hashing ids collapse to their base name. + + Hashing clones are left out. A collapsed base that is not already on the + panel is appended, so edgelist counts under that name are kept. + """ + hashing = {str(marker) for marker in hashing_antibodies} + kept = [str(marker) for marker in panel.markers if str(marker) not in hashing] + present = set(kept) + for marker in sorted(hashing): + base = collapsed_hashing_marker_id(marker) + if base not in present: + kept.append(base) + present.add(base) + return kept + + def _build_post_sample_calling_anndata( con: duckdb.DuckDBPyConnection, old_adata: anndata.AnnData, @@ -235,14 +253,15 @@ def _build_post_sample_calling_anndata( old_adata, hashing_antibody_mapping.hashing_antibodies ) - # Create the anndata object and remove all panel hashing markers from var - new_adata = pna_edgelist_to_anndata(con, panel) - non_hashing_markers = [ - marker - for marker in new_adata.var.index - if marker not in hashing_antibody_mapping.hashing_antibodies - ] - new_adata = new_adata[:, non_hashing_markers].copy() + # The edgelist already uses collapsed hashing ids (B2M-1 -> B2M). Include + # a collapsed base that the panel does not define, and leave the clones out. + new_adata = pna_edgelist_to_anndata( + con, + panel, + markers=_count_markers_after_hash_collapse( + panel, hashing_antibody_mapping.hashing_antibodies + ), + ) new_adata = add_missing_adata_info(new_adata, old_adata) # `sample` is a reserved, transient column added when reading a @@ -345,7 +364,7 @@ def sample_calling( enrichment_threshold, undetermined_sample_name, ) - panel = PNAAntibodyPanel.from_pxl_dataset(input_pxl) + panel = read_dataset_panel(input_pxl) dehashed = hash_info.group_by("called_sample") output_files: list[Path] = [] @@ -360,8 +379,7 @@ def sample_calling( "sample_name": sample_name, "version": __version__, "technology": "single-cell-pna", - "panel_name": panel.name, - "panel_version": panel.version, + "hashing_collapsed": True, } nodes_to_remove = _find_nodes_to_remove( @@ -390,7 +408,7 @@ def sample_calling( { "hashed_marker": hashed_markers, "base_marker": [ - re.sub(r"-\d+$", "", marker) for marker in hashed_markers + collapsed_hashing_marker_id(marker) for marker in hashed_markers ], } ) @@ -435,6 +453,7 @@ def sample_calling( with PixelFileWriter(target_path) as pxl_file_writer: pxl_file_writer.write_metadata(metadata) + pxl_file_writer.write_panel(panel) pxl_file_writer.write_edgelist(Path(tmp_edgelist_parquet.name)) pxl_file_writer.write_adata(adata) output_files.append(target_path) diff --git a/tests/common/data_generator/molecules.py b/tests/common/data_generator/molecules.py index e58d175c0..7b56c0567 100644 --- a/tests/common/data_generator/molecules.py +++ b/tests/common/data_generator/molecules.py @@ -11,12 +11,20 @@ import numpy as np import polars as pl +from pixelator.pna.config.panel import split_hashing_marker_id from tests.common.data_generator.topology import generate_cell_graph if TYPE_CHECKING: from pixelator.pna.config.panel import PNAAntibodyPanel +def _hashing_index(marker_id: str) -> int: + parts = split_hashing_marker_id(marker_id) + if parts is None: + raise ValueError(f"Hashing marker {marker_id!r} must end with -.") + return int(parts[1]) + + def generate_edgelist( n_cells: int, n_nodes: int, @@ -197,7 +205,7 @@ def _hashing_indices_per_cell( if hashing.size == 0: return [None] * n_cells if hashing_indices is None: - indices = np.unique([int(m.rsplit("-", 1)[-1]) for m in hashing]) + indices = np.unique([_hashing_index(str(m)) for m in hashing]) else: indices = np.unique(hashing_indices) cell_indices = np.resize(indices, n_cells) @@ -324,7 +332,7 @@ def _assign_markers( return node_umi_map hashing_markers = markers[is_hashing] - index = np.array([int(m.rsplit("-", 1)[-1]) for m in hashing_markers]) + index = np.array([_hashing_index(str(m)) for m in hashing_markers]) chosen = hashing_markers[index == hashing_index] if chosen.size == 0: return node_umi_map diff --git a/tests/common/data_generator/pixelfile.py b/tests/common/data_generator/pixelfile.py index 520c5fb8c..41546a597 100644 --- a/tests/common/data_generator/pixelfile.py +++ b/tests/common/data_generator/pixelfile.py @@ -55,12 +55,11 @@ def write_pna_pxl( "sample_name": sample_name, "version": __version__, "technology": "single-cell-pna", - "panel_name": panel.name, - "panel_version": panel.version, } ) writer.write_edgelist(pxl_edgelist) adata = pna_edgelist_to_anndata(writer.get_connection(), panel) + writer.write_panel(panel) writer.write_adata(adata) return path diff --git a/tests/common/data_generator/test_molecules.py b/tests/common/data_generator/test_molecules.py index 003ffa261..fe0814a79 100644 --- a/tests/common/data_generator/test_molecules.py +++ b/tests/common/data_generator/test_molecules.py @@ -11,6 +11,7 @@ from polars.testing import assert_frame_equal from pixelator.pna.config import pna_config +from pixelator.pna.config.panel import split_hashing_marker_id from tests.common.data_generator.molecules import ( _assign_markers, _assign_umis, @@ -154,7 +155,7 @@ def test_assign_markers_hashing_index(marker_panel): by_index = defaultdict(list) for marker, hashing in zip(markers, is_hashing): if hashing: - by_index[int(marker.rsplit("-", 1)[-1])].append(marker) + by_index[int(split_hashing_marker_id(marker)[1])].append(marker) # the chosen index's markers share the extra ~f and sit above the base chosen = by_index[hashing_index] diff --git a/tests/common/data_generator/test_pixelfile.py b/tests/common/data_generator/test_pixelfile.py index e386cf6a6..a55c1d624 100644 --- a/tests/common/data_generator/test_pixelfile.py +++ b/tests/common/data_generator/test_pixelfile.py @@ -9,6 +9,7 @@ from pixelator.pna.config import pna_config from pixelator.pna.config.panel import load_antibody_panel from pixelator.pna.pixeldataset import read +from pixelator.pna.pixeldataset.io import read_dataset_panel from tests.common.data_generator import ( generate_edgelist, write_pna_pxl, @@ -102,11 +103,12 @@ def test_write_pna_pxl_adata_metrics(written_pxl): def test_write_pna_pxl_metadata(written_pxl, real_panel): - """Metadata records the sample name, technology and panel identity.""" - metadata = read(written_pxl).metadata() - - assert set(metadata) == {SAMPLE_NAME} - entry = metadata[SAMPLE_NAME] + """Sample metadata does not duplicate the panel tables.""" + dataset = read(written_pxl) + entry = dataset.metadata()[SAMPLE_NAME] assert entry["technology"] == "single-cell-pna" - assert entry["panel_name"] == real_panel.name - assert entry["panel_version"] == real_panel.version + assert "panel_name" not in entry + assert "panel_version" not in entry + panel = read_dataset_panel(dataset) + assert panel.name == real_panel.name + assert panel.version == real_panel.version diff --git a/tests/pna/config/test_config.py b/tests/pna/config/test_config.py index 458131c41..b925379c6 100644 --- a/tests/pna/config/test_config.py +++ b/tests/pna/config/test_config.py @@ -9,6 +9,7 @@ import pandas as pd import pytest +from pixelator.common.config import AntibodyPanelMetadata from pixelator.pna.config.assay import PNARegionType, get_position_in_parent from pixelator.pna.config.config_class import ( PNAConfig, @@ -19,6 +20,31 @@ from pixelator.pna.config.panel import PNAAntibodyPanel, load_antibody_panel +def test_add_panel_rejects_multiple_sources(): + """A config panel is one CSV, so a concatenation cannot be registered.""" + + def panel(name: str, marker: str, sequence: str) -> PNAAntibodyPanel: + frame = pd.DataFrame( + [ + { + "marker_id": marker, + "control": False, + "sequence_1": sequence, + "sequence_2": sequence, + } + ] + ).set_index("marker_id") + return PNAAntibodyPanel( + frame, AntibodyPanelMetadata(name=name, version="1.0.0") + ) + + combined = PNAAntibodyPanel.concatenate( + [panel("base", "CD3", "AAAA"), panel("addon", "CD19", "CCCC")] + ) + with pytest.raises(ValueError, match="single-source"): + PNAConfig().add_panel(combined) + + def test_config_creation(): """Verify config creation.""" config = PNAConfig() diff --git a/tests/pna/config/test_panel.py b/tests/pna/config/test_panel.py index 6a2b7416b..6a4de66ed 100644 --- a/tests/pna/config/test_panel.py +++ b/tests/pna/config/test_panel.py @@ -11,6 +11,11 @@ from pixelator.common.config import AntibodyPanelMetadata from pixelator.pna.config.panel import PNAAntibodyPanel from pixelator.pna.pixeldataset import read +from pixelator.pna.pixeldataset.io import read_dataset_panel + + +def _one_metadata() -> AntibodyPanelMetadata: + return AntibodyPanelMetadata(name="test_panel", version="0.0.0") @pytest.fixture @@ -65,7 +70,7 @@ def test_panel_properties(panel_df): Args: panel_df: panel df. """ - panel = PNAAntibodyPanel(df=panel_df, metadata=None) + panel = PNAAntibodyPanel(df=panel_df, metadata=_one_metadata()) def test_panel_validation_fails_on_underscores_in_marker_names(panel_df): @@ -80,7 +85,7 @@ def test_panel_validation_fails_on_underscores_in_marker_names(panel_df): AssertionError, match=r".*The marker_id column should not contain underscores.*Offending values:.*", ): - PNAAntibodyPanel(df=panel_df, metadata=None) + PNAAntibodyPanel(df=panel_df, metadata=_one_metadata()) def test_panel_validation_fails_on_white_space_in_marker_names(panel_df): @@ -95,7 +100,7 @@ def test_panel_validation_fails_on_white_space_in_marker_names(panel_df): AssertionError, match=r".*The marker_id column should not contain white-spaces.*Offending values:.*", ): - PNAAntibodyPanel(df=panel_df, metadata=None) + PNAAntibodyPanel(df=panel_df, metadata=_one_metadata()) def test_panel_validation_fails_on_invalid_uniprot_ids(panel_df): @@ -110,7 +115,7 @@ def test_panel_validation_fails_on_invalid_uniprot_ids(panel_df): AssertionError, match=r".*Invalid UniProt IDs found.*Please conform to the naming convention or remove the following IDs:.*", ): - PNAAntibodyPanel(df=panel_df, metadata=None) + PNAAntibodyPanel(df=panel_df, metadata=_one_metadata()) def test_panel_validation_ok_on_concatenated_uniprot_ids(panel_df): @@ -120,7 +125,7 @@ def test_panel_validation_ok_on_concatenated_uniprot_ids(panel_df): panel_df: panel df. """ panel_df.loc["marker1", "uniprot_id"] = "P05107;P15391" - PNAAntibodyPanel(df=panel_df, metadata=None) + PNAAntibodyPanel(df=panel_df, metadata=_one_metadata()) def test_panel_validation_ok_uniprotid_empty(panel_df): @@ -130,7 +135,7 @@ def test_panel_validation_ok_uniprotid_empty(panel_df): panel_df: panel df. """ panel_df.loc["marker1", "uniprot_id"] = "" - PNAAntibodyPanel(df=panel_df, metadata=None) + PNAAntibodyPanel(df=panel_df, metadata=_one_metadata()) def test_panel_from_pxl(pxl_file): @@ -139,7 +144,7 @@ def test_panel_from_pxl(pxl_file): Args: pxl_file: pxl file. """ - panel = PNAAntibodyPanel.from_pxl_dataset(read(pxl_file)) + panel = read_dataset_panel(read(pxl_file)) assert panel.name == "test-pna-panel" assert panel.version == "0.1.0" assert panel.description == "Test R&D panel for RNA" diff --git a/tests/pna/config/test_panel_concatenation.py b/tests/pna/config/test_panel_concatenation.py new file mode 100644 index 000000000..a2ad8bcfc --- /dev/null +++ b/tests/pna/config/test_panel_concatenation.py @@ -0,0 +1,347 @@ +"""Tests for concatenating several antibody panels into one. + +Copyright © 2026 Pixelgen Technologies AB. +""" + +from pathlib import Path + +import numpy as np +import pandas as pd +import pytest +from anndata import AnnData +from pandas.testing import assert_frame_equal + +from pixelator.common.config import AntibodyPanelMetadata +from pixelator.pna.config.panel import ( + PNAAntibodyPanel, + align_panel_patches, + aligned_dataset_panel, +) +from pixelator.pna.pixeldataset import read +from pixelator.pna.pixeldataset.io import PixelFileWriter, PxlFile, read_dataset_panel + + +def _panel( + name: str, version: str, rows: list[dict], product: str = "proxiome" +) -> PNAAntibodyPanel: + frame = pd.DataFrame(rows).set_index("marker_id") + return PNAAntibodyPanel( + frame, + AntibodyPanelMetadata(name=name, version=version, product=product), + file_name=f"{name}.csv", + ) + + +def _marker(marker_id: str, sequence: str, **extra) -> dict: + row = { + "marker_id": marker_id, + "control": False, + "sequence_1": sequence, + "sequence_2": sequence, + } + row.update(extra) + return row + + +def test_concatenate_keeps_metadata_for_each_source(): + base = _panel( + "base", + "1.0.0", + [_marker("CD3", "AAAA")], + ) + addon = _panel( + "addon", + "1.2.0", + [_marker("CD19", "CCCC")], + ) + + combined = PNAAntibodyPanel.concatenate([base, addon]) + + for field in ( + "metadata", + "name", + "version", + "product", + "description", + "aliases", + "archived", + ): + with pytest.raises(ValueError, match="single source"): + getattr(combined, field) + assert [source.metadata.name for source in combined.sources] == ["base", "addon"] + assert [source.metadata.version for source in combined.sources] == [ + "1.0.0", + "1.2.0", + ] + assert list(combined.markers) == ["CD3", "CD19"] + assert combined.marker_source_ids.loc["CD3"] == 0 + assert combined.marker_source_ids.loc["CD19"] == 1 + + +def test_equality_ignores_row_order_and_filename_but_not_marker_source(): + left = _panel("base", "1.0.0", [_marker("CD3", "AAAA"), _marker("CD19", "CCCC")]) + right = PNAAntibodyPanel( + left.df.iloc[::-1], + AntibodyPanelMetadata(name="base", version="1.0.0", product="proxiome"), + file_name="other.csv", + ) + assert left == right + + base = _panel("base", "1.0.0", [_marker("CD3", "AAAA")]) + addon = _panel("addon", "1.2.0", [_marker("CD19", "CCCC")]) + frame = pd.concat([base.df, addon.df]) + frame.index.name = "marker_id" + sources = [*base.sources, *addon.sources] + assigned = PNAAntibodyPanel( + frame, + sources=sources, + marker_source_ids=pd.Series({"CD3": 0, "CD19": 1}, dtype="int64"), + ) + swapped = PNAAntibodyPanel( + frame, + sources=sources, + marker_source_ids=pd.Series({"CD3": 1, "CD19": 0}, dtype="int64"), + ) + assert assigned != swapped + + reordered = PNAAntibodyPanel.concatenate([addon, base]) + assert PNAAntibodyPanel.concatenate([base, addon]) == reordered + + +def test_aligned_dataset_panel_ignores_source_order(): + base_old = _panel("base", "1.0.0", [_marker("CD3", "AAAA")], product="kit") + base_new = _panel("base", "1.0.1", [_marker("CD3E", "AAAA")], product="kit") + addon = _panel("addon", "2.0.0", [_marker("CD19", "CCCC")], product="kit") + forward = PNAAntibodyPanel.concatenate([base_old, addon]) + reverse = PNAAntibodyPanel.concatenate([addon, base_new]) + + assert aligned_dataset_panel([forward, reverse]) == PNAAntibodyPanel.concatenate( + [base_new, addon] + ) + + +def test_concatenate_one_panel_keeps_its_metadata(): + panel = _panel("base", "1.0.0", [_marker("CD3", "AAAA")]) + assert PNAAntibodyPanel.concatenate([panel]) is panel + assert panel.name == "base" + + +def test_concatenate_blanks_optional_columns_missing_from_one_source(): + with_uniprot = _panel( + "base", + "1.0.0", + [_marker("CD3", "AAAA", uniprot_id="P01730")], + ) + without_uniprot = _panel("addon", "1.0.0", [_marker("CD19", "CCCC")]) + + combined = PNAAntibodyPanel.concatenate([with_uniprot, without_uniprot]) + + assert combined.df.loc["CD3", "uniprot_id"] == "P01730" + assert combined.df.loc["CD19", "uniprot_id"] == "" + + hashing = _panel( + "hashing", + "1.0.0", + [_marker("B2M-1", "GGGG", sample_hashing=True)], + ) + combined_hashing = PNAAntibodyPanel.concatenate([without_uniprot, hashing]) + assert combined_hashing.hashing_marker_ids == {"B2M-1"} + + +def test_replace_source_blanks_optional_columns_missing_from_one_side(): + old_base = _panel("base", "1.0.0", [_marker("CD3", "AAAA")], product="kit") + new_base = _panel( + "base", + "1.0.1", + [_marker("CD3E", "AAAA", uniprot_id="P07766")], + product="kit", + ) + addon = _panel("addon", "2.0.0", [_marker("CD19", "CCCC")], product="kit") + old = PNAAntibodyPanel.concatenate([old_base, addon]) + new = PNAAntibodyPanel.concatenate([new_base, addon]) + + replaced = old.replace_source(0, new_base) + + assert replaced.df.loc["CD3E", "uniprot_id"] == "P07766" + assert replaced.df.loc["CD19", "uniprot_id"] == "" + assert replaced == new + assert aligned_dataset_panel([old, new]) == new + + addon_with_uniprot = _panel( + "addon", + "2.0.0", + [_marker("CD19", "CCCC", uniprot_id="P15391")], + product="kit", + ) + newer_without = _panel("base", "1.0.1", [_marker("CD3E", "AAAA")], product="kit") + kept_column = PNAAntibodyPanel.concatenate([old_base, addon_with_uniprot]) + replaced_other_side = kept_column.replace_source(0, newer_without) + assert replaced_other_side.df.loc["CD3E", "uniprot_id"] == "" + assert replaced_other_side.df.loc["CD19", "uniprot_id"] == "P15391" + + +def test_patch_bump_matches_when_the_added_column_is_not_last(): + old_base = _panel("base", "1.0.0", [_marker("CD3", "AAAA")], product="kit") + new_frame = pd.DataFrame([_marker("CD3E", "AAAA")]).set_index("marker_id") + new_frame["uniprot_id"] = "P07766" + new_frame = new_frame[["control", "uniprot_id", "sequence_1", "sequence_2"]] + new_base = PNAAntibodyPanel( + new_frame, + AntibodyPanelMetadata(name="base", version="1.0.1", product="kit"), + file_name="base.csv", + ) + addon = _panel("addon", "2.0.0", [_marker("CD19", "CCCC")], product="kit") + old = PNAAntibodyPanel.concatenate([old_base, addon]) + new = PNAAntibodyPanel.concatenate([new_base, addon]) + updated, _, _ = align_panel_patches([old, new]) + + assert list(updated[0].df.columns) != list(updated[1].df.columns) + assert updated[0] == updated[1] + assert aligned_dataset_panel([old, new]) == new + + +def test_concatenate_rejects_duplicate_marker_and_sequence(): + left = _panel("base", "1.0.0", [_marker("CD3", "AAAA")]) + same_marker = _panel("addon", "1.0.0", [_marker("CD3", "CCCC")]) + with pytest.raises(AssertionError, match="marker_id were not unique"): + PNAAntibodyPanel.concatenate([left, same_marker]) + + same_sequence = _panel("addon", "1.0.0", [_marker("CD19", "AAAA")]) + with pytest.raises(AssertionError, match="sequence_1"): + PNAAntibodyPanel.concatenate([left, same_sequence]) + + +def test_hashing_marker_ids_need_a_numeric_suffix_and_must_not_nest(): + missing_suffix = pd.DataFrame( + [_marker("B2M", "AAAA", sample_hashing=True)] + ).set_index("marker_id") + with pytest.raises(AssertionError, match="must end with -"): + PNAAntibodyPanel( + missing_suffix, + AntibodyPanelMetadata(name="base", version="1.0.0"), + ) + + nested = pd.DataFrame( + [ + _marker("B2M-1", "AAAA", sample_hashing=True), + _marker("B2M-1-1", "CCCC", sample_hashing=True), + ] + ).set_index("marker_id") + with pytest.raises(AssertionError, match="must not collapse"): + PNAAntibodyPanel( + nested, + AntibodyPanelMetadata(name="base", version="1.0.0"), + ) + + +def test_panel_tables_roundtrip_and_legacy_uns(tmp_path: Path): + panel = _panel( + "base", + "1.0.0", + [_marker("CD3", "AAAA"), _marker("CD19", "CCCC", control=True)], + ) + path = tmp_path / "sample.pxl" + with PixelFileWriter(path) as writer: + writer.write_metadata({"sample_name": "sample"}) + writer.write_panel(panel) + # A minimal AnnData so the file is still a pxl for readers that only + # need the panel tables. + connection = writer.get_connection() + connection.execute("CREATE TABLE edgelist (umi1 INTEGER)") + connection.execute('CREATE TABLE "__adata__X" (index VARCHAR)') + connection.execute('CREATE TABLE "__adata__obs" (index VARCHAR)') + connection.execute('CREATE TABLE "__adata__var" (index VARCHAR)') + + loaded = PxlFile(path).read_panel() + assert loaded is not None + assert loaded.name == "base" + assert loaded.version == "1.0.0" + assert loaded.sources[0].metadata.product == "proxiome" + assert_frame_equal(loaded.df, panel.df) + + legacy = tmp_path / "legacy.pxl" + adata = AnnData(X=np.zeros((1, panel.df.shape[0])), var=panel.df.copy()) + adata.uns["panel_metadata"] = panel.metadata.model_dump() + adata.uns["panel_metadata"]["panel_columns"] = list(panel.df.columns) + with PixelFileWriter(legacy) as writer: + writer.write_metadata({"sample_name": "legacy"}) + writer.write_adata(adata) + legacy_panel = PxlFile(legacy).read_panel() + assert legacy_panel is not None + assert legacy_panel.name == "base" + assert_frame_equal(legacy_panel.df, panel.df) + + +def test_patch_bump_is_per_source_and_skips_collapsed_hashing_clones(): + old_base = _panel("base", "1.0.0", [_marker("CD3", "AAAA")], product="kit") + new_base = _panel("base", "1.0.1", [_marker("CD3E", "AAAA")], product="kit") + addon = _panel("addon", "2.0.0", [_marker("CD19", "CCCC")], product="kit") + old = PNAAntibodyPanel.concatenate([old_base, addon]) + new = PNAAntibodyPanel.concatenate([new_base, addon]) + + old_adata = AnnData(var=pd.DataFrame(index=["CD3", "CD19"])) + new_adata = AnnData(var=pd.DataFrame(index=["CD3E", "CD19"])) + upgraded, renames, _hash_renames = align_panel_patches( + [old, new], [old_adata, new_adata] + ) + + assert upgraded[0].sources[0].metadata.version == "1.0.1" + assert upgraded[0].sources[1].metadata.version == "2.0.0" + assert "CD3E" in upgraded[0].markers + assert renames[0] == {"CD3": "CD3E"} + assert renames[1] == {} + + hashing = _panel( + "hash", + "1.0.0", + [ + _marker("B2M-1", "GGGG", sample_hashing=True), + _marker("CD3", "AAAA"), + ], + product="kit", + ) + hashing_next = _panel( + "hash", + "1.0.1", + [ + _marker("B2M-1", "GGGG", sample_hashing=True), + _marker("CD3E", "AAAA"), + ], + product="kit", + ) + collapsed = AnnData(var=pd.DataFrame(index=["B2M", "CD3"])) + current = AnnData(var=pd.DataFrame(index=["B2M-1", "CD3E"])) + with pytest.raises(ValueError, match="Missing markers"): + align_panel_patches( + [hashing, hashing_next], + [collapsed, current], + pxl_file_metadata=[{"hashing_collapsed": False}, {}], + ) + upgraded, renames, _hash_renames = align_panel_patches( + [hashing, hashing_next], + [collapsed, current], + pxl_file_metadata=[{"hashing_collapsed": True}, {}], + ) + assert renames[0]["CD3"] == "CD3E" + assert "B2M-1" in upgraded[0].markers + + +def test_read_dataset_panel_rejects_a_file_from_before_0_22(tmp_path: Path): + path = tmp_path / "old.pxl" + with PixelFileWriter(path) as writer: + writer.write_metadata({"sample_name": "old", "panel_name": "base"}) + connection = writer.get_connection() + connection.execute("CREATE TABLE edgelist (umi1 INTEGER)") + connection.execute('CREATE TABLE "__adata__X" (index VARCHAR)') + connection.execute('CREATE TABLE "__adata__obs" (index VARCHAR)') + connection.execute('CREATE TABLE "__adata__var" (index VARCHAR)') + connection.execute('CREATE TABLE "__adata__uns" (value JSON)') + + with pytest.raises(ValueError, match="current version of the software"): + read_dataset_panel(read(path)) + + +def test_pxl_fixture_roundtrip(pxl_file): + panel = read_dataset_panel(read(pxl_file)) + assert panel.name == "test-pna-panel" + assert "panel_metadata" not in read(pxl_file).adata().uns diff --git a/tests/pna/conftest.py b/tests/pna/conftest.py index 9043bef63..a313d08fd 100644 --- a/tests/pna/conftest.py +++ b/tests/pna/conftest.py @@ -297,12 +297,12 @@ def create_pxl_file( writer.write_edgelist(edgelist_parquet_path) con = writer.get_connection() adata = pna_edgelist_to_anndata(con, panel=panel) + writer.write_panel(panel) writer.write_adata(adata) writer.write_metadata( { "sample_name": sample_name, "version": "0.1.0", - "panel_name": "custom_panel", } ) if proximity_parquet_path: diff --git a/tests/pna/denoise/test_denoise.py b/tests/pna/denoise/test_denoise.py index f345d61ff..9cbc283bd 100644 --- a/tests/pna/denoise/test_denoise.py +++ b/tests/pna/denoise/test_denoise.py @@ -4,7 +4,6 @@ """ from pathlib import Path -from unittest import mock import networkx as nx import numpy as np @@ -213,22 +212,10 @@ def test_denoise_one_core_analysis(synthetic_denoise_pxl_dataset, tmp_path): pxl_file_target = PixelDatasetSaver(pxl_dataset=synthetic_denoise_pxl_dataset).save( "PNA055_Sample07_S7", Path(tmp_path) / "layout.pxl" ) - with mock.patch( - "pixelator.pna.analysis.denoise.load_antibody_panel" - ) as mock_load_panel: - # This is a workaround to make sure that the correct panel is loaded - # eventhough we no longer set a default panel file. - def f(*args, **kwargs): - return load_antibody_panel(pna_config, "proxiome-v1-immuno-155-v1.0") - - mock_load_panel.side_effect = f - - manager = AnalysisManager( - [DenoiseGraph(run_one_core=True, run_ace=False)], n_cores=1 - ) - denoised_dataset = manager.execute( - synthetic_denoise_pxl_dataset, pxl_file_target - ) + manager = AnalysisManager( + [DenoiseGraph(run_one_core=True, run_ace=False)], n_cores=1 + ) + denoised_dataset = manager.execute(synthetic_denoise_pxl_dataset, pxl_file_target) adata = denoised_dataset.adata() obs = adata.obs @@ -308,16 +295,10 @@ def _run_one_core_denoise(dataset, target_path): target = PixelDatasetSaver(pxl_dataset=dataset).save( "PNA055_Sample07_S7", target_path ) - with mock.patch( - "pixelator.pna.analysis.denoise.load_antibody_panel" - ) as mock_load_panel: - mock_load_panel.side_effect = lambda *args, **kwargs: load_antibody_panel( - pna_config, "proxiome-v1-immuno-155-v1.0" - ) - manager = AnalysisManager( - [DenoiseGraph(run_one_core=True, run_ace=False)], n_cores=1 - ) - return manager.execute(dataset, target) + manager = AnalysisManager( + [DenoiseGraph(run_one_core=True, run_ace=False)], n_cores=1 + ) + return manager.execute(dataset, target) def test_denoise_one_core_analysis_deterministic( @@ -503,18 +484,10 @@ def test_denoise_ace_analysis(synthetic_denoise_pxl_dataset, tmp_path): "PNA055_Sample07_S7", Path(tmp_path) / "layout.pxl" ) - with mock.patch( - "pixelator.pna.analysis.denoise.load_antibody_panel" - ) as mock_load_panel: - mock_load_panel.side_effect = lambda *args, **kwargs: load_antibody_panel( - pna_config, "proxiome-v1-immuno-155-v1.0" - ) - manager = AnalysisManager( - [DenoiseGraph(run_one_core=False, run_ace=True)], n_cores=1 - ) - denoised_dataset = manager.execute( - synthetic_denoise_pxl_dataset, pxl_file_target - ) + manager = AnalysisManager( + [DenoiseGraph(run_one_core=False, run_ace=True)], n_cores=1 + ) + denoised_dataset = manager.execute(synthetic_denoise_pxl_dataset, pxl_file_target) obs = denoised_dataset.adata().obs @@ -557,18 +530,10 @@ def test_denoise_ace_pls_one_core(synthetic_denoise_pxl_dataset, tmp_path): "PNA055_Sample07_S7", Path(tmp_path) / "layout.pxl" ) - with mock.patch( - "pixelator.pna.analysis.denoise.load_antibody_panel" - ) as mock_load_panel: - mock_load_panel.side_effect = lambda *args, **kwargs: load_antibody_panel( - pna_config, "proxiome-v1-immuno-155-v1.0" - ) - manager = AnalysisManager( - [DenoiseGraph(run_one_core=True, run_ace=True, run_pls=True)], n_cores=1 - ) - denoised_dataset = manager.execute( - synthetic_denoise_pxl_dataset, pxl_file_target - ) + manager = AnalysisManager( + [DenoiseGraph(run_one_core=True, run_ace=True, run_pls=True)], n_cores=1 + ) + denoised_dataset = manager.execute(synthetic_denoise_pxl_dataset, pxl_file_target) obs = denoised_dataset.adata().obs diff --git a/tests/pna/pixeldataset/test_anndata_helper.py b/tests/pna/pixeldataset/test_anndata_helper.py index e5daf6056..b5c64747c 100644 --- a/tests/pna/pixeldataset/test_anndata_helper.py +++ b/tests/pna/pixeldataset/test_anndata_helper.py @@ -3,15 +3,18 @@ Copyright © 2025 Pixelgen Technologies AB. """ +import shutil from pathlib import Path +import duckdb import numpy as np import polars as pl import pytest from pixelator.common.utils.testing import adata_assert_equal -from pixelator.pna.config.panel import PNAAntibodyPanel +from pixelator.pna.config.panel import PNAAntibodyPanel, align_panel_patches from pixelator.pna.pixeldataset import PNAPixelDataset +from pixelator.pna.pixeldataset.io import Query, read_dataset_panel from pixelator.pna.pixeldataset.io.anndata_helper import AnnDataHelper from tests.pna.conftest import create_pxl_file @@ -70,6 +73,8 @@ def _build_two_sample_dataset_with_panels( edgelist_parquet_path: Path, panel_old: PNAAntibodyPanel, panel_new: PNAAntibodyPanel, + proximity_old: Path | None = None, + proximity_new: Path | None = None, ) -> PNAPixelDataset: """Create two on-disk PXL samples with distinct panels for bumping patch version tests. @@ -83,7 +88,7 @@ def _build_two_sample_dataset_with_panels( target=tmp_path / "sample_old.pxl", sample_name="sample_old", edgelist_parquet_path=edgelist_parquet_path, - proximity_parquet_path=None, + proximity_parquet_path=proximity_old, layout_parquet_path=None, panel=panel_old, ) @@ -100,7 +105,7 @@ def _build_two_sample_dataset_with_panels( target=tmp_path / "sample_new.pxl", sample_name="sample_new", edgelist_parquet_path=sample_new_edgelist, - proximity_parquet_path=None, + proximity_parquet_path=proximity_new, layout_parquet_path=None, panel=panel_new, ) @@ -111,16 +116,20 @@ class TestAnnDataHelper: """Represent test ann data helper.""" def test_anndata_helper_matches_dataset_adata_no_transforms( - self, pxl_dataset, adata_data + self, pxl_dataset, adata_data, panel ): """Verify anndata helper matches dataset adata no transforms. Args: pxl_dataset: pxl dataset. adata_data: adata data. + panel: panel. """ adata_data = adata_data.copy() adata_data.obs["sample"] = "test_sample" + adata_data.var = adata_data.var.join( + panel.df.reindex(adata_data.var_names), how="left" + ) helper = AnnDataHelper(pxl_dataset.view) res = helper.read_adata(add_clr_transform=False, add_log1p_transform=False) @@ -245,20 +254,32 @@ def test_bumps_to_latest_patch_when_prerequisites_are_met( positive_cells_count = np.random.randint(0, 100, adata_old.var.shape[0]) adata_old.var["positive_cells_count"] = positive_cells_count - assert adata_old.var.loc["MarkerA", "uniprot_id"] == "P12345" - assert adata_new.var.loc["MarkerANew", "uniprot_id"] == "Q9UPN0" + assert "MarkerA" in adata_old.var_names + assert "MarkerANew" in adata_new.var_names + assert "uniprot_id" not in adata_old.var.columns assert "target_class" not in adata_old.var.columns - assert "target_class" in adata_new.var.columns + adata_old.var["uniprot_id"] = "STALE" + adata_old.var["retired_field"] = "old-value" + adata_old.uns["panel_metadata"] = { + "name": "stale", + "panel_columns": ["uniprot_id", "retired_field"], + } - bumped = helper._try_bump_adata_panel_version([adata_old, adata_new]) + bumped = dataset.view.apply_panel_patch_to_adatas([adata_old, adata_new]) - assert "target_class" in bumped[0].var.columns - assert bumped[0].var.loc["MarkerANew", "target_class"] == "new-value" + assert "MarkerANew" in bumped[0].var_names + assert "MarkerA" not in bumped[0].var_names assert bumped[0].var.loc["MarkerANew", "uniprot_id"] == "Q9UPN0" - - assert "target_class" in bumped[1].var.columns + assert bumped[0].var.loc["MarkerANew", "target_class"] == "new-value" + assert "retired_field" not in bumped[0].var.columns assert bumped[1].var.loc["MarkerANew", "uniprot_id"] == "Q9UPN0" assert bumped[1].var.loc["MarkerANew", "target_class"] == "new-value" + assert "panel_metadata" not in bumped[0].uns + assert "panel_metadata" not in bumped[1].uns + + aligned = read_dataset_panel(dataset) + assert aligned.df.loc["MarkerANew", "uniprot_id"] == "Q9UPN0" + assert aligned.df.loc["MarkerANew", "target_class"] == "new-value" assert "positive_cells_count" in bumped[0].var.columns assert "positive_cells_count" not in bumped[1].var.columns @@ -270,6 +291,251 @@ def test_bumps_to_latest_patch_when_prerequisites_are_met( assert (adata_old[:, "MarkerC"].X == bumped[0][:, "MarkerC"].X).all() assert (adata_new[:, "MarkerC"].X == bumped[1][:, "MarkerC"].X).all() + def test_edgelist_view_exposes_bumped_marker_ids( + self, + tmp_path: Path, + edgelist_parquet_path: Path, + panel: PNAAntibodyPanel, + ): + """The session edgelist uses bumped ids before any later query.""" + panel_old = _panel_with_version_product_and_uniprot( + panel, + version="0.1.0", + product="test-product", + marker_a_uniprot="P12345", + ) + panel_new = _panel_with_version_product_and_uniprot( + panel, + version="0.1.1", + product="test-product", + marker_a_uniprot="Q9UPN0", + marker_a_new_name="MarkerANew", + ) + dataset = _build_two_sample_dataset_with_panels( + tmp_path=tmp_path, + edgelist_parquet_path=edgelist_parquet_path, + panel_old=panel_old, + panel_new=panel_new, + ) + with dataset.view.open() as session: + frame = session.execute_eager( + Query("SELECT sample, marker_1, marker_2 FROM edgelist", {}) + ) + old = frame.filter(pl.col("sample") == "sample_old") + old_ids = set(old["marker_1"].to_list() + old["marker_2"].to_list()) + assert "MarkerA" not in old_ids + assert "MarkerANew" in old_ids + + def test_sample_filter_keeps_the_patch_bump( + self, + tmp_path: Path, + edgelist_parquet_path: Path, + panel: PNAAntibodyPanel, + ): + """A cohort narrowed to one sample still uses the bumped marker ids.""" + panel_old = _panel_with_version_product_and_uniprot( + panel, + version="0.1.0", + product="test-product", + marker_a_uniprot="P12345", + ) + panel_new = _panel_with_version_product_and_uniprot( + panel, + version="0.1.1", + product="test-product", + marker_a_uniprot="Q9UPN0", + marker_a_new_name="MarkerANew", + ) + dataset = _build_two_sample_dataset_with_panels( + tmp_path=tmp_path, + edgelist_parquet_path=edgelist_parquet_path, + panel_old=panel_old, + panel_new=panel_new, + ) + + filtered = dataset.filter(samples={"sample_old"}) + adata = filtered.adata(add_clr_transform=False, add_log1p_transform=False) + assert "MarkerANew" in adata.var_names + assert "MarkerA" not in adata.var_names + edgelist = filtered.edgelist().to_polars() + ids = set(edgelist["marker_1"].to_list() + edgelist["marker_2"].to_list()) + assert "MarkerA" not in ids + assert "MarkerANew" in ids + + narrowed = dataset.filter(samples={"sample_old"}, markers={"MarkerANew"}) + narrowed_adata = narrowed.adata( + add_clr_transform=False, add_log1p_transform=False + ) + assert list(narrowed_adata.var_names) == ["MarkerANew"] + assert narrowed_adata.n_obs > 0 + + alone = PNAPixelDataset.from_pxl_files([tmp_path / "sample_old.pxl"]) + alone_adata = alone.adata(add_clr_transform=False, add_log1p_transform=False) + assert "MarkerA" in alone_adata.var_names + assert "MarkerANew" not in alone_adata.var_names + + def test_failed_patch_bump_leaves_the_cache_unset( + self, + tmp_path: Path, + edgelist_parquet_path: Path, + panel: PNAAntibodyPanel, + monkeypatch: pytest.MonkeyPatch, + ): + """A failed bump can be read again and then uses the new marker ids.""" + panel_old = _panel_with_version_product_and_uniprot( + panel, + version="0.1.0", + product="test-product", + marker_a_uniprot="P12345", + ) + panel_new = _panel_with_version_product_and_uniprot( + panel, + version="0.1.1", + product="test-product", + marker_a_uniprot="Q9UPN0", + marker_a_new_name="MarkerANew", + ) + dataset = _build_two_sample_dataset_with_panels( + tmp_path=tmp_path, + edgelist_parquet_path=edgelist_parquet_path, + panel_old=panel_old, + panel_new=panel_new, + ) + attempts = {"count": 0} + + def fail_once(*args, **kwargs): + attempts["count"] += 1 + if attempts["count"] == 1: + raise ValueError("bump failed") + return align_panel_patches(*args, **kwargs) + + monkeypatch.setattr( + "pixelator.pna.pixeldataset.io.pixel_data_viewer.align_panel_patches", + fail_once, + ) + with pytest.raises(ValueError, match="bump failed"): + dataset.view.marker_renames_by_sample() + assert dataset.view._marker_renames is None + assert dataset.view._hash_count_renames is None + assert dataset.view._upgraded_panels is None + + adata = dataset.adata(add_clr_transform=False, add_log1p_transform=False) + assert "MarkerANew" in adata.var_names + assert "MarkerA" not in adata.var_names + + def test_proximity_filter_uses_renamed_marker_ids( + self, + tmp_path: Path, + edgelist_parquet_path: Path, + panel: PNAAntibodyPanel, + ): + """A filter on the bumped name still finds rows stored under the old id.""" + panel_old = _panel_with_version_product_and_uniprot( + panel, + version="0.1.0", + product="test-product", + marker_a_uniprot="P12345", + ) + panel_new = _panel_with_version_product_and_uniprot( + panel, + version="0.1.1", + product="test-product", + marker_a_uniprot="Q9UPN0", + marker_a_new_name="MarkerANew", + ) + old_proximity = tmp_path / "old_proximity.parquet" + new_proximity = tmp_path / "new_proximity.parquet" + pl.DataFrame( + { + "component": ["fc07dea9b679aca7", "fc07dea9b679aca7"], + "marker_1": ["MarkerA", "MarkerB"], + "marker_2": ["MarkerA", "MarkerC"], + } + ).write_parquet(old_proximity) + pl.DataFrame( + { + "component": ["fc07dea9b679aca7_sample_new"], + "marker_1": ["MarkerANew"], + "marker_2": ["MarkerANew"], + } + ).write_parquet(new_proximity) + dataset = _build_two_sample_dataset_with_panels( + tmp_path=tmp_path, + edgelist_parquet_path=edgelist_parquet_path, + panel_old=panel_old, + panel_new=panel_new, + proximity_old=old_proximity, + proximity_new=new_proximity, + ) + + proximity = dataset.filter(markers={"MarkerANew"}).proximity( + add_marker_counts=False, add_logratio=False + ) + frame = proximity.to_polars() + pairs = set( + zip(frame["marker_1"].to_list(), frame["marker_2"].to_list(), strict=True) + ) + assert pairs == {("MarkerANew", "MarkerANew")} + assert len(proximity) == 2 + + def test_record_batches_use_renamed_marker_ids( + self, + tmp_path: Path, + edgelist_parquet_path: Path, + panel: PNAAntibodyPanel, + ): + """The streamed edgelist uses the same marker ids as ``to_polars``.""" + panel_old = _panel_with_version_product_and_uniprot( + panel, + version="0.1.0", + product="test-product", + marker_a_uniprot="P12345", + ) + panel_new = _panel_with_version_product_and_uniprot( + panel, + version="0.1.1", + product="test-product", + marker_a_uniprot="Q9UPN0", + marker_a_new_name="MarkerANew", + ) + dataset = _build_two_sample_dataset_with_panels( + tmp_path=tmp_path, + edgelist_parquet_path=edgelist_parquet_path, + panel_old=panel_old, + panel_new=panel_new, + ) + streamed = pl.concat( + [pl.from_arrow(batch) for batch in dataset.edgelist().to_record_batches()], + how="vertical", + ) + loaded = dataset.edgelist().to_polars() + assert streamed.sort(streamed.columns).equals(loaded.sort(loaded.columns)) + assert ( + "MarkerANew" + in streamed["marker_1"].to_list() + streamed["marker_2"].to_list() + ) + + def test_record_batches_rename_legacy_marker_columns( + self, tmp_path: Path, pxl_file: Path + ): + """A legacy edgelist stream uses the same column names as ``to_polars``.""" + legacy = tmp_path / "legacy.pxl" + shutil.copy(pxl_file, legacy) + with duckdb.connect(str(legacy)) as connection: + connection.execute("ALTER TABLE edgelist RENAME COLUMN marker_1 TO marker1") + connection.execute("ALTER TABLE edgelist RENAME COLUMN marker_2 TO marker2") + dataset = PNAPixelDataset.from_files(legacy) + streamed = pl.concat( + [pl.from_arrow(batch) for batch in dataset.edgelist().to_record_batches()], + how="vertical", + ) + loaded = dataset.edgelist().to_polars() + assert "marker_1" in streamed.columns + assert "marker_2" in streamed.columns + assert "marker1" not in streamed.columns + assert "marker2" not in streamed.columns + assert streamed.sort(streamed.columns).equals(loaded.sort(loaded.columns)) + @pytest.mark.parametrize( "new_version,new_product", [ @@ -325,9 +591,14 @@ def test_skips_bump_when_prerequisites_are_not_met( session=session, sample="sample_new" ) - not_bumped = helper._try_bump_adata_panel_version([adata_old, adata_new]) + not_bumped = dataset.view.apply_panel_patch_to_adatas([adata_old, adata_new]) + assert "MarkerA" in not_bumped[0].var_names + assert not_bumped[0].var.loc["MarkerA", "uniprot_id"] == "P12345" assert "target_class" not in not_bumped[0].var.columns + assert not_bumped[1].var.loc["MarkerA", "uniprot_id"] == "Q9UPN0" + assert not_bumped[1].var.loc["MarkerA", "target_class"] == "new-version" + assert "panel_metadata" not in not_bumped[0].uns assert (adata_old[:, "MarkerC"].X == not_bumped[0][:, "MarkerC"].X).all() assert (adata_new[:, "MarkerC"].X == not_bumped[1][:, "MarkerC"].X).all() diff --git a/tests/pna/pixeldataset/test_pixel_file_writer.py b/tests/pna/pixeldataset/test_pixel_file_writer.py index 8b3d1a390..85f93be0b 100644 --- a/tests/pna/pixeldataset/test_pixel_file_writer.py +++ b/tests/pna/pixeldataset/test_pixel_file_writer.py @@ -1,11 +1,17 @@ """Copyright © 2025 Pixelgen Technologies AB.""" +import json from pathlib import Path +import numpy as np +import pandas as pd import polars as pl import pytest +from anndata import AnnData -from pixelator.pna.pixeldataset.io import PixelFileWriter +from pixelator.common.config import AntibodyPanelMetadata +from pixelator.pna.config.panel import PNAAntibodyPanel +from pixelator.pna.pixeldataset.io import PixelFileWriter, PxlFile class TestPixelFileWriter: @@ -129,3 +135,127 @@ def test_write_layouts(self, tmp_path, layout_parquet_path): target = tmp_path / "file.pxl" with PixelFileWriter(target) as writer: writer.write_layouts(layout_parquet_path) + + def test_write_adata_omits_panel_columns_already_stored_in_the_panel_tables( + self, tmp_path + ): + """Panel fields joined onto var in memory stay out of the written var.""" + panel = PNAAntibodyPanel( + pd.DataFrame( + { + "control": [False], + "sequence_1": ["AAAA"], + "sequence_2": ["AAAA"], + "uniprot_id": ["P01730"], + }, + index=pd.Index(["CD3"], name="marker_id"), + ), + AntibodyPanelMetadata(name="base", version="1.0.0", product="proxiome"), + ) + adata = AnnData( + X=np.zeros((1, 1)), + obs=pd.DataFrame(index=["c1"]), + var=pd.DataFrame( + { + "antibody_count": [1], + "control": [False], + "sequence_1": ["AAAA"], + "sequence_2": ["AAAA"], + "uniprot_id": ["P01730"], + }, + index=pd.Index(["CD3"], name="marker_id"), + ), + ) + target = tmp_path / "file.pxl" + with PixelFileWriter(target) as writer: + writer.write_panel(panel) + writer.write_adata(adata) + columns = set( + writer.get_connection() + .execute("SELECT * FROM __adata__var LIMIT 0") + .fetchdf() + .columns + ) + assert "antibody_count" in columns + assert "control" not in columns + assert "sequence_1" not in columns + assert "sequence_2" not in columns + assert "uniprot_id" not in columns + + def test_write_adata_keeps_var_columns_when_no_panel_is_stored(self, tmp_path): + """A file without panel tables keeps whatever columns var already has.""" + adata = AnnData( + X=np.zeros((1, 1)), + obs=pd.DataFrame(index=["c1"]), + var=pd.DataFrame( + {"antibody_count": [1], "control": [False]}, + index=pd.Index(["CD3"], name="marker_id"), + ), + ) + target = tmp_path / "file.pxl" + with PixelFileWriter(target) as writer: + writer.write_adata(adata) + columns = set( + writer.get_connection() + .execute("SELECT * FROM __adata__var LIMIT 0") + .fetchdf() + .columns + ) + assert "antibody_count" in columns + assert "control" in columns + + def test_write_adata_stores_a_legacy_panel_before_replacing_it(self, tmp_path): + """A 0.22–0.30 file keeps its panel when AnnData is written again.""" + legacy = AnnData( + X=np.zeros((1, 1)), + obs=pd.DataFrame(index=["c1"]), + var=pd.DataFrame( + { + "antibody_count": [1], + "control": [False], + "sequence_1": ["AAAA"], + "sequence_2": ["AAAA"], + }, + index=pd.Index(["CD3"], name="marker_id"), + ), + ) + legacy.uns["panel_metadata"] = { + "name": "base", + "version": "1.0.0", + "product": "proxiome", + "panel_columns": ["control", "sequence_1", "sequence_2"], + } + target = tmp_path / "file.pxl" + with PixelFileWriter(target) as writer: + writer.write_adata(legacy) + + rewritten = AnnData( + X=np.zeros((1, 1)), + obs=pd.DataFrame({"k_core_1": [3]}, index=["c1"]), + var=legacy.var.copy(), + ) + with PixelFileWriter(target) as writer: + writer.write_adata(rewritten) + columns = set( + writer.get_connection() + .execute("SELECT * FROM __adata__var LIMIT 0") + .fetchdf() + .columns + ) + uns_row = ( + writer.get_connection() + .execute("SELECT value FROM __adata__uns") + .fetchone() + ) + + uns = json.loads(uns_row[0]) if isinstance(uns_row[0], str) else uns_row[0] + panel = PxlFile(target).read_panel() + assert panel is not None + assert panel.name == "base" + assert panel.version == "1.0.0" + assert list(panel.markers) == ["CD3"] + assert "panel_metadata" not in uns + assert "antibody_count" in columns + assert "control" not in columns + assert "sequence_1" not in columns + assert "sequence_2" not in columns diff --git a/tests/pna/pixeldataset/test_pna_pixeldataset.py b/tests/pna/pixeldataset/test_pna_pixeldataset.py index b23e87eb4..b2eed6071 100644 --- a/tests/pna/pixeldataset/test_pna_pixeldataset.py +++ b/tests/pna/pixeldataset/test_pna_pixeldataset.py @@ -129,15 +129,19 @@ def test_view(self, pxl_dataset): """ assert isinstance(pxl_dataset.view, PixelDataViewer) - def test_adata(self, pxl_dataset: PNAPixelDataset, adata_data): + def test_adata(self, pxl_dataset: PNAPixelDataset, adata_data, panel): """Verify adata. Args: adata_data: adata data. pxl_dataset: Pxl dataset. + panel: panel. """ adata_data = adata_data.copy() adata_data.obs["sample"] = "test_sample" + adata_data.var = adata_data.var.join( + panel.df.reindex(adata_data.var_names), how="left" + ) adata = pxl_dataset.adata(add_clr_transform=False, add_log1p_transform=False) adata_assert_equal(adata, adata_data) @@ -187,7 +191,7 @@ def test_metadata_returns_expected_mapping(self, pxl_dataset: PNAPixelDataset): assert metadata.keys() == {"test_sample"} assert metadata["test_sample"]["sample_name"] == "test_sample" assert metadata["test_sample"]["version"] == "0.1.0" - assert metadata["test_sample"]["panel_name"] == "custom_panel" + assert "panel_name" not in metadata["test_sample"] def test_metadata_returns_empty_dict_when_metadata_table_is_empty( self, diff --git a/tests/pna/sample_calling/test_sample_calling.py b/tests/pna/sample_calling/test_sample_calling.py index f8752eba5..45a28c1e0 100644 --- a/tests/pna/sample_calling/test_sample_calling.py +++ b/tests/pna/sample_calling/test_sample_calling.py @@ -15,7 +15,7 @@ from pixelator.common.config import AntibodyPanelMetadata from pixelator.pna import read from pixelator.pna.anndata import pna_edgelist_to_anndata -from pixelator.pna.config.panel import PNAAntibodyPanel +from pixelator.pna.config.panel import PNAAntibodyPanel, collapsed_hashing_marker_id from pixelator.pna.pixeldataset.io import PixelFileWriter from pixelator.pna.sample_calling import ( _add_original_hash_counts_to_obs, @@ -385,8 +385,10 @@ def test_sample_calling(sample_hashed_pixel_files, tmp_path): dehashed_counts = dehashed_adata.to_df() assert dehashed_counts.shape[0] == 10 sample_name = list(dehashed_pxl.metadata().keys())[0] + assert dehashed_pxl.metadata()[sample_name]["hashing_collapsed"] is True + assert "sample_calling" not in (dehashed_adata.uns or {}) for ab in hashed_antibodies[sample_name]: - base_name = ab.split("-")[0] + base_name = collapsed_hashing_marker_id(ab) assert all( original_counts.loc[ dehashed_counts.index, @@ -484,12 +486,12 @@ def test_sample_calling_does_not_strip_suffix_from_non_hash_markers( writer.write_edgelist(edgelist) con = writer.get_connection() adata = pna_edgelist_to_anndata(con, panel=panel) + writer.write_panel(panel) writer.write_adata(adata) writer.write_metadata( { "sample_name": "input", "version": "0.1.0", - "panel_name": "custom_panel", } ) @@ -520,6 +522,71 @@ def test_sample_calling_does_not_strip_suffix_from_non_hash_markers( assert dehashed_counts.loc["c1", "PD-1"] == original_counts.loc["c1", "PD-1"] +def test_sample_calling_keeps_collapsed_hash_marker_missing_from_panel( + tmp_path: Path, +): + """A collapsed hashing id stays in var when the panel has no such marker.""" + panel_df = pd.DataFrame( + [ + { + "marker_id": "PD-1", + "control": False, + "uniprot_id": "P00001", + "sequence_1": "ATCGATCGAA", + "sequence_2": "ATCGATCGAC", + }, + { + "marker_id": "HashA-1", + "control": False, + "uniprot_id": "P00003", + "sequence_1": "ATCGATCGTT", + "sequence_2": "ATCGATCGTG", + }, + ] + ).set_index("marker_id") + panel = PNAAntibodyPanel( + df=panel_df, + metadata=AntibodyPanelMetadata(name="test-panel", version="0.1.0"), + ) + edgelist = pl.DataFrame( + { + "umi1": pl.Series([1, 2, 3], dtype=pl.UInt64), + "umi2": pl.Series([2, 3, 4], dtype=pl.UInt64), + "read_count": pl.Series([10, 10, 10], dtype=pl.UInt32), + "uei_count": pl.Series([5, 5, 5], dtype=pl.UInt32), + "marker_1": ["PD-1", "PD-1", "HashA-1"], + "marker_2": ["PD-1", "HashA-1", "HashA-1"], + "component": ["c1", "c1", "c1"], + } + ) + target = tmp_path / "input.pxl" + with PixelFileWriter(target) as writer: + writer.write_edgelist(edgelist) + con = writer.get_connection() + adata = pna_edgelist_to_anndata(con, panel=panel) + writer.write_panel(panel) + writer.write_adata(adata) + writer.write_metadata({"sample_name": "input", "version": "0.1.0"}) + + out_dir = tmp_path / "out" + out_dir.mkdir() + sample_calling( + input_pxl=read(target), + hashing_antibody_mapping=HashedAntibodyMapping( + mapping={"S1": ["HashA-1"]}, + all_hashing_antibodies=["HashA-1"], + ), + output_folder=out_dir, + remove_incompatible=True, + enrichment_threshold=1.5, + ) + + dehashed_counts = read(next(out_dir.glob("*.dehashed.pxl"))).adata().to_df() + assert "HashA" not in panel.markers + assert "HashA-1" not in dehashed_counts.columns + assert dehashed_counts.loc["c1", "HashA"] == 3 + + @pytest.mark.slow def test_sample_calling_without_uei_count(tmp_path: Path): """sample_calling succeeds when the input edgelist has no uei_count column.""" @@ -583,12 +650,12 @@ def test_sample_calling_without_uei_count(tmp_path: Path): writer.write_edgelist(edgelist) con = writer.get_connection() adata = pna_edgelist_to_anndata(con, panel=panel) + writer.write_panel(panel) writer.write_adata(adata) writer.write_metadata( { "sample_name": "input", "version": "0.1.0", - "panel_name": "custom_panel", } ) @@ -760,7 +827,9 @@ def adata( n = len(self._hash_enrichment_factors) # Explicit string obs index avoids AnnData ImplicitModificationWarning on index coercion. obs_index = [f"comp_{i}" for i in range(n)] - obs = {"hash_enrichment_factor": self._hash_enrichment_factors} + obs: dict[str, list[float] | list[int]] = { + "hash_enrichment_factor": self._hash_enrichment_factors + } if self._reads_in_component is not None: assert len(self._reads_in_component) == n obs["reads_in_component"] = self._reads_in_component diff --git a/tests/pna/test_anndata.py b/tests/pna/test_anndata.py index 2f3369a20..3edc69a50 100644 --- a/tests/pna/test_anndata.py +++ b/tests/pna/test_anndata.py @@ -155,8 +155,6 @@ def create_pixel_dataset_connection(edgelist): "sample_name": "temp", "version": __version__, "technology": "single-cell-pna", - "panel_name": "mock.name", - "panel_version": "mock.version", } ) writer.write_edgelist(edgelist.collect()) @@ -192,28 +190,16 @@ def test_pna_edgelist_to_anndata(pixelconnection): "antibody_count": 7, "antibody_pct": np.float32(0.4375), "components": 3, - "uniprot_id": "P61769", - "control": False, - "nuclear": "yes", - "sequence_1": "AAAA", }, "B": { "antibody_count": 5, "antibody_pct": np.float32(0.3125), "components": 3, - "uniprot_id": "P05107", - "control": False, - "nuclear": "no", - "sequence_1": "CCCC", }, "C": { "antibody_count": 4, "antibody_pct": np.float32(0.25), "components": 2, - "uniprot_id": "", - "control": True, - "nuclear": "no", - "sequence_1": "GGGG", }, }, orient="index",