diff --git a/scripts/create_resources/combine/process_datasets_xenium_nebius.sh b/scripts/create_resources/combine/process_datasets_xenium_nebius.sh index a630053f1..259208fd5 100644 --- a/scripts/create_resources/combine/process_datasets_xenium_nebius.sh +++ b/scripts/create_resources/combine/process_datasets_xenium_nebius.sh @@ -129,16 +129,6 @@ param_list: dataset_description: "Xenium V1 FFPE Human Breast IDC + 2021 Wu scRNAseq" dataset_organism: "homo_sapiens" - - id: "2026_10x_human_breast_cancer_atera_combined" - input_sp: "$input_dir/10x_atera/2026_10x_human_breast_cancer_atera/dataset.zarr" - input_sc: "$input_dir/wu_human_breast_cancer_sc/2021Wu_human_breast_cancer_sc/dataset.h5ad" - dataset_id: "2026_10x_human_breast_cancer_atera_combined" - dataset_name: "Human breast cancer combined 2026 10x Atera WTA 2021 Wu scRNAseq" - dataset_url: "https://www.10xgenomics.com/datasets/atera-wta-ffpe-human-breast-cancer" - dataset_reference: "https://doi.org/10.1038/s41588-021-00911-1" - dataset_summary: "Atera WTA FFPE Human Breast Cancer (DCIS Grade 3) + 2021 Wu scRNAseq" - dataset_description: "Atera WTA FFPE Human Breast Cancer (DCIS Grade 3) + 2021 Wu scRNAseq" - dataset_organism: "homo_sapiens" output_sc: "\$id/output_sc.h5ad" output_sp: "\$id/output_sp.zarr" diff --git a/scripts/run_benchmark/cellposev4_params.yaml b/scripts/run_benchmark/cellposev4_params.yaml new file mode 100644 index 000000000..6681b0f89 --- /dev/null +++ b/scripts/run_benchmark/cellposev4_params.yaml @@ -0,0 +1,47 @@ +# Parameter sweep for the cellposev4 (Cellpose 4 / Cellpose-SAM) segmentation method. +# Shared source of truth for both run_test_cellposev4_local.sh (read as a local file) +# and run_test_cellposev4_nebius.sh (read from GitHub via a raw URL, since the Nebius +# compute env pulls the repo but cannot see the launch host's local files). +# +# Consumed by the run_benchmark workflow via the `method_parameters_yaml` setting +# (src/workflows/run_benchmark/main.nf). For every method the workflow builds: +# * one "default" variant using the `default:` args below, and +# * one extra variant per value in each `sweep:` list, with that ONE arg overridden. +# The benchmark varies a SINGLE parameter at a time (a "star" around the default, not +# a full grid), so total cellposev4 variants = 1 default + sum(sweep list lengths) = 19. +# +# See src/methods_segmentation/cellposev4/NOTES.md ("Optimization / tuning") for the +# rationale, tiers, and defaults. Trim or widen the lists to taste. +parameters: + cellposev4: + # Baseline == the component's shipped (deliberately speed-tuned) defaults. + default: + diameter: 30.0 + flow_threshold: 0.0 + cellprob_threshold: 0.0 + niter: 10 + min_size: -1 + resample: false + sweep: + # ---- Tier 1: highest impact on segmentation quality ---- + # diameter: cellpose rescales so objects land near ~30 px (training mean), + # so 30.0 == "no rescaling". Xenium morphology (~0.2125 um/px) => nuclei + # ~35-40 px, whole cells ~60-70 px. Span small nuclei -> whole cells. + diameter: [15.0, 20.0, 40.0, 60.0] + # cellprob_threshold: recall<->precision dial (network logit ~ -6..+6, def 0). + # Lower -> recover more/dimmer cells; higher -> drop dim-region detections. + cellprob_threshold: [-3.0, -1.0, 1.0, 3.0] + # flow_threshold: 0.0 (our default) disables flow-shape QC (keeps ill-shaped + # masks); >0 re-enables it, more permissive as it rises. 0.4 = Cellpose default. + flow_threshold: [0.4, 0.6, 0.8] + # ---- Tier 2: quality/speed trade-offs ---- + # niter: flow-dynamics iterations; 10 (default) is aggressively low and may + # under-converge large/elongated cells. Cellpose's own default is ~diameter- + # proportional (often ~200). + niter: [50, 100, 200] + # resample: BOOLEAN, default false -> the only meaningful non-default value is + # true (run dynamics at full resolution: smoother boundaries, slower). + resample: [true] + # min_size: -1 (default) disables small-mask removal; positive values drop + # progressively larger specks/debris. + min_size: [15, 30, 50] diff --git a/src/datasets/loaders/bruker_cosmx/config.vsh.yaml b/src/datasets/loaders/bruker_cosmx/config.vsh.yaml index f51200615..901efe036 100644 --- a/src/datasets/loaders/bruker_cosmx/config.vsh.yaml +++ b/src/datasets/loaders/bruker_cosmx/config.vsh.yaml @@ -4,10 +4,15 @@ namespace: datasets/loaders argument_groups: - name: Inputs arguments: - - type: file + - type: string name: --input_raw example: "https://smi-public.objects.liquidweb.services/HalfBrain.zip" - description: "Download file url for the raw data" + description: | + URL/path to the raw data zip. Passed as a string (not a staged file) on + purpose: the script streams it in place and extracts only the ~15 GB sopa + needs (Morphology2D + CellLabels), skipping the ~190 GB of Morphology3D / + AnalysisResults. Staging the full zip (175 GiB) alongside the extraction + overruns the node's scratch and OOMs. Supports s3:// (via s3fs) or a local path. - type: file name: --input_flat_files example: "https://smi-public.objects.liquidweb.services/Half%20%20Brain%20simple%20%20files%20.zip" @@ -65,9 +70,10 @@ engines: __merge__: - /src/base/setup_spatialdata_partial.yaml setup: - - type: python + - type: python pypi: - sopa + - s3fs # stream the raw zip from s3:// and extract only what sopa needs (see --input_raw) - type: native runners: diff --git a/src/datasets/loaders/bruker_cosmx/script.py b/src/datasets/loaders/bruker_cosmx/script.py index 00b2e3f64..883913f29 100644 --- a/src/datasets/loaders/bruker_cosmx/script.py +++ b/src/datasets/loaders/bruker_cosmx/script.py @@ -199,6 +199,59 @@ def extract_zip(input_zip: Path, output_dir: Path, strip_root: bool = False): else: raise AssertionError(f"No CellLabels folder and no per-FOV CellLabels tifs found in {DATA_DIR}") +############################################# +# Memory optimization: lean transcript read # +############################################# +# sopa's CosMx reader loads the ENTIRE `*_tx_file.csv` into a pandas DataFrame in +# one `pd.read_csv` (sopa/io/reader/cosmx.py::_CosMXReader.read_transcripts). For +# the mouse-brain half-hemisphere that's hundreds of millions of rows, and the +# string columns ('target', 'CellComp', 'cell', ...) come in as object dtype +# (~60 B/row) — this is the peak that blows past even the 480 GB veryhighmem cap. +# +# We wrap the `read_csv` used inside sopa's cosmx module so that, *only for the +# transcripts file*, we (a) read the 'target' feature column as categorical and +# (b) keep just the columns sopa's stitched-FOV path (read_transcripts + +# _get_global_cell_id) and the downstream transcript schema actually use. Every +# other sopa read (fov positions, metadata, counts, polygons) is left untouched, +# and all of sopa's coordinate/stitching math runs unchanged on the leaner frame. +from sopa.io.reader import cosmx as _cosmx_mod + +_real_pd = _cosmx_mod.pd +_real_read_csv = _real_pd.read_csv + +# Columns used by read_transcripts (fov=None) + _get_global_cell_id, plus 'z' +# which the downstream transcript schema allows. Everything else in the tx file +# (x_local_px, y_local_px, cell, CellComp, ...) is unused downstream. +_TX_KEEP = ["fov", "cell_ID", "target", "x_global_px", "y_global_px", "z"] +_TX_REQUIRED = {"target", "x_global_px", "y_global_px", "fov", "cell_ID"} + +def _lean_read_csv(filepath_or_buffer, *args, **kwargs): + name = str(filepath_or_buffer) + if "_tx_file.csv" in name and "usecols" not in kwargs: + header = _real_read_csv( + filepath_or_buffer, nrows=0, + compression=kwargs.get("compression", "infer"), + ) + cols = set(header.columns) + if _TX_REQUIRED.issubset(cols): + kwargs["usecols"] = [c for c in _TX_KEEP if c in cols] + dtype = dict(kwargs.get("dtype") or {}) + dtype.setdefault("target", "category") + kwargs["dtype"] = dtype + log(f"Lean transcript read of {Path(name).name}: usecols={kwargs['usecols']}, target=category") + return _real_read_csv(filepath_or_buffer, *args, **kwargs) + +class _PandasReadCsvProxy: + """Forwards every attribute to pandas, but leans the transcripts read_csv.""" + def __init__(self, real): + self.__dict__["_real"] = real + def read_csv(self, *args, **kwargs): + return _lean_read_csv(*args, **kwargs) + def __getattr__(self, item): + return getattr(self._real, item) + +_cosmx_mod.pd = _PandasReadCsvProxy(_real_pd) + ######################################### # Convert raw files to spatialdata zarr # ######################################### diff --git a/src/datasets/workflows/process_bruker_cosmx/config.vsh.yaml b/src/datasets/workflows/process_bruker_cosmx/config.vsh.yaml index 3c256c603..195b96788 100644 --- a/src/datasets/workflows/process_bruker_cosmx/config.vsh.yaml +++ b/src/datasets/workflows/process_bruker_cosmx/config.vsh.yaml @@ -4,10 +4,10 @@ namespace: datasets/workflows argument_groups: - name: Inputs arguments: - - type: file + - type: string name: --input_raw example: "https://smi-public.objects.liquidweb.services/HalfBrain.zip" - description: "Download file url for the raw data" + description: "URL/path to the raw data zip. A string (not a staged file) so the loader can stream it and extract only what sopa needs — see the loader's --input_raw." - type: file name: --input_flat_files example: "https://smi-public.objects.liquidweb.services/Half%20%20Brain%20simple%20%20files%20.zip" diff --git a/src/methods_cell_type_annotation/rctd/script.R b/src/methods_cell_type_annotation/rctd/script.R index 371f1556e..1bdc3e6c4 100644 --- a/src/methods_cell_type_annotation/rctd/script.R +++ b/src/methods_cell_type_annotation/rctd/script.R @@ -18,7 +18,7 @@ par <- list( ) meta <- list( - 'cpus': 4, + cpus = 4 ) ## VIASH END diff --git a/src/methods_segmentation/cellposev4/config.vsh.yaml b/src/methods_segmentation/cellposev4/config.vsh.yaml index 25d20d264..7b994cf7d 100644 --- a/src/methods_segmentation/cellposev4/config.vsh.yaml +++ b/src/methods_segmentation/cellposev4/config.vsh.yaml @@ -12,7 +12,11 @@ links: documentation: "https://cellpose.readthedocs.io/en/latest/" repository: "https://github.com/mouseland/cellpose" references: - doi: "10.1038/s41592-020-01018-x" + doi: + # Cellpose-SAM (cpsam): the model this component actually runs + - "10.1101/2025.04.28.651001" + # Original Cellpose: the flow-dynamics framework it builds on + - "10.1038/s41592-020-01018-x" __merge__: /src/api/comp_method_segmentation.yaml @@ -25,6 +29,10 @@ arguments: type: double default: 0.0 description: "Flow error threshold. Set to 0 to skip flow quality check for faster execution." + - name: --cellprob_threshold + type: double + default: 0.0 + description: "Cell probability threshold (network output ~-6 to 6). Pixels above it are used to seed masks. Lower to recover more/dimmer cells; raise to drop detections in dim/low-contrast regions." - name: --niter type: integer default: 10 diff --git a/src/methods_segmentation/cellposev4/script.py b/src/methods_segmentation/cellposev4/script.py index 1c5e51a7f..159584419 100644 --- a/src/methods_segmentation/cellposev4/script.py +++ b/src/methods_segmentation/cellposev4/script.py @@ -18,6 +18,7 @@ "output": "segmentation.zarr", "diameter": 30.0, "flow_threshold": 0.0, + "cellprob_threshold": 0.0, "niter": 10, "min_size": -1, "resample": False, @@ -50,7 +51,7 @@ def convert_to_lower_dtype(arr): print('Initializing Cellpose model', flush=True) model = CellposeModel(gpu=torch.cuda.is_available()) -eval_params = {k: par[k] for k in ("diameter", "flow_threshold", "niter", "min_size", "resample") if par.get(k) is not None} +eval_params = {k: par[k] for k in ("diameter", "flow_threshold", "cellprob_threshold", "niter", "min_size", "resample") if par.get(k) is not None} print('Running Cellpose segmentation with parameters:', eval_params, flush=True) masks, _, _ = model.eval(image[0], progress=True, **eval_params) diff --git a/src/methods_transcript_assignment/comseg/script.py b/src/methods_transcript_assignment/comseg/script.py index 1905fe29b..84ed9a438 100644 --- a/src/methods_transcript_assignment/comseg/script.py +++ b/src/methods_transcript_assignment/comseg/script.py @@ -37,22 +37,37 @@ # Read input files print('Reading input files', flush=True) -sdata = sd.read_zarr(par['input_ist']) +sdata_src = sd.read_zarr(par['input_ist']) sdata_segm = sd.read_zarr(par['input_segmentation']) # Multi-partition parquet restarts its index at 0 per partition, producing duplicate # index labels that break sopa's cell_id assignment and the final write with # "cannot reindex on an axis with duplicate labels". Rebuild the transcripts as a # single-partition frame with a clean RangeIndex before any sopa op touches them. -_tx = sdata[par['transcripts_key']] +_tx = sdata_src[par['transcripts_key']] _tx_reset = dd.from_pandas(_tx.compute().reset_index(drop=True), npartitions=1) _tx_reset.attrs.update(_tx.attrs) -del sdata[par['transcripts_key']] -sdata[par['transcripts_key']] = _tx_reset # Convert the prior segmentation to polygons -sdata["segmentation_boundaries"] = sd.to_polygons(sdata_segm["segmentation"]) -del sdata["segmentation_boundaries"]["label"] # make_transcript_patches will create a new label column and fails if one exists. +segmentation_boundaries = sd.to_polygons(sdata_segm["segmentation"]) +del segmentation_boundaries["label"] # make_transcript_patches will create a new label column and fails if one exists. + +# Run sopa/ComSeg on a FRESH, in-memory SpatialData. sd.read_zarr() leaves `sdata_src` +# backed by the shared source zarr (sdata_src.path -> input_ist), and sopa's +# make_image_patches / make_transcript_patches / segmentation.comseg all write their new +# elements (image_patches, transcripts_patches, comseg_boundaries) and a .sopa_cache back +# through that path -- mutating the supposed-to-be-immutable source dataset in place. If a +# concurrent run or an OOM kill interrupts one of those writes, it leaves a half-written +# element (e.g. shapes/image_patches with no shapes.parquet) that then breaks EVERY +# downstream reader of the dataset. baysor and proseg avoid this by running sopa on a fresh +# in-memory object; do the same here (carrying the image for make_image_patches and the +# prior boundaries for the transcript patches) so all sopa writes land in the ephemeral +# work dir, never the source zarr. +sdata = sd.SpatialData( + images={"image": sdata_src["image"]}, + points={par['transcripts_key']: _tx_reset}, + shapes={"segmentation_boundaries": segmentation_boundaries}, +) # Make patches sopa.make_image_patches(sdata, image_key="image", patch_width=par["patch_width"], patch_overlap=par["patch_overlap"]) diff --git a/src/methods_transcript_assignment/segger/config.vsh.yaml b/src/methods_transcript_assignment/segger/config.vsh.yaml index 889e143b9..f83d61b0e 100644 --- a/src/methods_transcript_assignment/segger/config.vsh.yaml +++ b/src/methods_transcript_assignment/segger/config.vsh.yaml @@ -116,14 +116,25 @@ engines: # supported numba/numpy range is a latent risk for segger's phenograph step. # Validated set: numba 0.60.0 / numpy 2.0.2 / llvmlite 0.43.0 (cudf .values, # cuml kNN, cugraph louvain, torch, spatialdata/anndata all import + run). + # - PIN polars>=1.24,<1.26 (== what cudf-polars 25.4 requires). segger's github + # install leaves polars unpinned; this keeps it at the RAPIDS-compatible 1.24. + # See the writer.py from_torch shim below for why we must stay <1.26. # - swap opencv-python -> opencv-python-headless: segger pulls opencv-python, # whose cv2 needs libGL.so.1 (a system lib we cannot apt-install as non-root); # the headless build provides the same cv2 without it. See NOTES.md. - type: docker run: - - pip install --no-cache-dir "pandas>=2.0,<2.2.4" "anndata>=0.12,<0.13" "scanpy>=1.11,<1.12" "spatialdata>=0.7,<0.8" "dask==2025.2.0" "distributed==2025.2.0" "numba>=0.59.1,<0.61" "numpy>=2.0,<2.1" + - pip install --no-cache-dir "pandas>=2.0,<2.2.4" "anndata>=0.12,<0.13" "scanpy>=1.11,<1.12" "spatialdata>=0.7,<0.8" "dask==2025.2.0" "distributed==2025.2.0" "numba>=0.59.1,<0.61" "numpy>=2.0,<2.1" "polars>=1.24,<1.26" - pip uninstall -y opencv-python || true - pip install --no-cache-dir opencv-python-headless + # segger's data/writer.py calls pl.from_torch to turn the prediction tensors into a + # polars frame, but polars only added from_torch AFTER 1.27 (PR pola-rs/polars#22177) + # -- and cudf-polars 25.4 caps polars at <1.26, so we cannot get the real one without + # breaking the RAPIDS stack. Inject the missing function (torch tensor -> from_numpy, + # honoring segger's list/dict schema) right after writer.py's `import torch`. getattr + # preserves a real from_torch if a future polars provides one. Path is fixed by the + # py3.12 base image; a build failure here (missing file) is a loud signal it moved. + - "sed -i '/^import torch$/a pl.from_torch = getattr(pl, \"from_torch\", lambda data, schema=None, **k: pl.from_numpy(data.detach().cpu().numpy(), schema=schema))' /opt/conda/lib/python3.12/site-packages/segger/data/writer.py" - type: native runners: