diff --git a/.dockerignore b/.dockerignore new file mode 100644 index 0000000..f841322 --- /dev/null +++ b/.dockerignore @@ -0,0 +1,18 @@ +# Keep the build context to what `uv sync --locked` actually needs. In +# particular .venv/ would otherwise ship a host-platform virtualenv into the +# image, and .claude/worktrees/ holds full copies of this repo. +# +# Do NOT add README.md or LICENSE here: pyproject.toml references them via +# `readme` and `license-files`, so the build fails without them. +.venv/ +.git/ +.claude/ +.pytest_cache/ +__pycache__/ +*.py[cod] +.coverage +coverage.xml +htmlcov/ +build/ +dist/ +*.egg-info/ diff --git a/CHANGELOG.md b/CHANGELOG.md index 1bd7f7d..2a129b5 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -3,6 +3,48 @@ Notable changes to `adata-cli`. Versions are `MAJOR.MINOR.PATCH`; tags carry no `v` prefix. +## 0.5.1 + +Makes the container image usable from Nextflow, and stops `copy_dataset` +reading one row at a time from row-chunked stores. + +### Fixed + +- **Copying a row-chunked store was dominated by read latency.** The read step + was the source's chunk height verbatim, so a store chunked `(1, n_cols)` was + copied one row per read. On a network filesystem (Lustre, NFS) each read is a + round-trip, so a million-row copy spent nearly all of its time waiting. Reads + are now sized to a 32 MiB budget, rounded down to a whole number of source + chunks. A `(1_000_000, 30_000)` float32 store chunked `(1, 30_000)` goes from + 1 row per read to 279. +- Read sizing no longer trusts `itemsize` for variable-length strings. h5py + reports 8 there because the value is a pointer, which overestimated the row + count by an order of magnitude and broke the memory bound. + +### Container + +- **The image could not be used from a Nextflow process.** Nextflow requires + `/bin/bash` to be the container entrypoint, so `ENTRYPOINT ["adata"]` made + every Docker- and Podman-backed task fail with `No such command + '/bin/bash'`. Apptainer was unaffected, as `singularity exec` ignores the + entrypoint. +- **Task metrics were silently lost.** `procps` is absent from the base image, + so Nextflow could not run `ps` to collect them. The required tool set + (`bash`, `ps`, `awk`, `date`, `grep`, `sed`, `tail`, `tee`) is now installed + and asserted at build time. +- `PYTHONNOUSERSITE` is set, so a bind-mounted `$HOME` under Apptainer can no + longer shadow the image's virtualenv with the user's `~/.local` packages. +- `XDG_CACHE_HOME` points at `/tmp`, so the image tolerates being run under an + arbitrary UID with no writable `$HOME`. +- Added a `.dockerignore`. Local builds were copying the host's `.venv`, + `.git` and `.pytest_cache` into the image. + +### Changed + +- **The image no longer sets an entrypoint, so the command must be named + explicitly:** `docker run IMAGE adata view file.h5ad`, where `docker run + IMAGE view file.h5ad` previously worked. + ## 0.5.0 Renamed from `h5ad` to `adata-cli`, restored compatibility with current diff --git a/Dockerfile b/Dockerfile index 5955f68..8d0e79b 100644 --- a/Dockerfile +++ b/Dockerfile @@ -1,24 +1,33 @@ # Base image: Python 3.12 + uv preinstalled (Debian slim) FROM ghcr.io/astral-sh/uv:python3.12-bookworm-slim -ENV UV_NO_DEV=1 +# PYTHONNOUSERSITE: Apptainer bind-mounts the host $HOME by default, so a user's +# ~/.local/lib/python3.12/site-packages would otherwise shadow this venv. +# XDG_CACHE_HOME: Nextflow is commonly configured with `-u $(id -u):$(id -g)`, +# which leaves the container with no writable $HOME. +# UV_COMPILE_BYTECODE: bake .pyc at build time, so nothing writes to a +# read-only rootfs on first import. +ENV UV_NO_DEV=1 \ + UV_COMPILE_BYTECODE=1 \ + PYTHONNOUSERSITE=1 \ + PYTHONDONTWRITEBYTECODE=1 \ + PYTHONUNBUFFERED=1 \ + XDG_CACHE_HOME=/tmp/.cache -WORKDIR /cli - -# Copy the project files (from the GitHub Actions checkout context) -COPY . . - -# --locked asserts that uv.lock is in sync with pyproject.toml, so an image -# can never be built from a lockfile that drifted. -RUN uv sync --locked +# procps supplies `ps`, which Nextflow needs to collect per-task metrics. +# curl, unzip and ca-certificates fetch duckdb below, and are left in place +# rather than purged: pipeline scripts routinely reach for curl, and TLS +# roots are worth having in any container that may touch the network. +RUN apt-get update \ + && apt-get install -y --no-install-recommends \ + ca-certificates curl procps mawk unzip \ + && rm -rf /var/lib/apt/lists/* # duckdb, for the filtering workflows in the docs: export obs to CSV, query it, # feed the names back to `adata subset --obs`. A single static binary, so it # needs no venv and cannot conflict with the project's dependencies. ARG DUCKDB_VERSION=v1.1.3 -RUN apt-get update \ - && apt-get install -y --no-install-recommends ca-certificates curl unzip \ - && ARCH="$(dpkg --print-architecture)" \ +RUN ARCH="$(dpkg --print-architecture)" \ && case "$ARCH" in \ amd64) DUCKDB_ARCH=amd64 ;; \ arm64) DUCKDB_ARCH=aarch64 ;; \ @@ -28,13 +37,43 @@ RUN apt-get update \ "https://github.com/duckdb/duckdb/releases/download/${DUCKDB_VERSION}/duckdb_cli-linux-${DUCKDB_ARCH}.zip" \ && unzip -q /tmp/duckdb.zip -d /usr/local/bin \ && chmod +x /usr/local/bin/duckdb \ - && rm /tmp/duckdb.zip \ - && apt-get purge -y curl unzip \ - && apt-get autoremove -y \ - && rm -rf /var/lib/apt/lists/* + && rm /tmp/duckdb.zip + +# Fail the build, rather than every Nextflow task, if the base image ever drops +# one of the tools Nextflow requires in a task container. +RUN set -eu; for t in bash ps awk date grep sed tail tee; do \ + command -v "$t" >/dev/null || { echo "missing required tool: $t" >&2; exit 1; }; \ + done + +WORKDIR /cli + +# Copy the project files (from the GitHub Actions checkout context) +COPY . . + +# --locked asserts that uv.lock is in sync with pyproject.toml, so an image +# can never be built from a lockfile that drifted. +# +# uv honours XDG_CACHE_HOME, so the sync leaves a root-owned package cache at +# /tmp/.cache -- which made the variable self-defeating, as a task running +# under an arbitrary UID then could not write to the very path it advertises. +# Clear it and leave an empty world-writable directory behind. Done in this +# same layer because a later `rm` would mask the files without reclaiming +# them; that reclaims about 9 MB, the cache being mostly hardlinks into the +# venv rather than separate copies. +RUN uv sync --locked \ + && rm -rf /tmp/.cache /tmp/uv-*.lock \ + && mkdir -p /tmp/.cache \ + && chmod 1777 /tmp/.cache # Put the project venv on PATH so `adata` is directly runnable ENV PATH="/cli/.venv/bin:${PATH}" -ENTRYPOINT ["adata"] -CMD ["--help"] +# No ENTRYPOINT on purpose: Nextflow requires /bin/bash to be the container +# entrypoint, so the image must not set one of its own. This is why invocations +# spell out the command: `docker run IMAGE adata view file.h5ad`. +# +# Deliberately NOT set here: OMP_NUM_THREADS / OPENBLAS_NUM_THREADS. NumPy's +# BLAS sizes its thread pool to the whole host, which oversubscribes a shared +# LSF node. This workload is streaming I/O, so capping it would cost nothing -- +# but it belongs in the pipeline's `env` scope, not baked into the image. +CMD ["adata", "--help"] diff --git a/README.md b/README.md index fb6c2af..1b55ad9 100644 --- a/README.md +++ b/README.md @@ -91,5 +91,5 @@ adata concat per_sample/*.h5ad -o merged.h5ad --join outer --label sample A docker image is available on QUAY: `quay.io/cellgeni/adata-cli:latest`. Pull and run with: ```bash -docker run --rm -it -v /path/to/data:/data quay.io/cellgeni/adata-cli:latest view /data/your_file.h5ad +docker run --rm -it -v /path/to/data:/data quay.io/cellgeni/adata-cli:latest adata view /data/your_file.h5ad ``` \ No newline at end of file diff --git a/docs/index.md b/docs/index.md index b37ac86..ed04b98 100644 --- a/docs/index.md +++ b/docs/index.md @@ -22,7 +22,7 @@ Or run it without installing anything: ```bash docker run --rm -it -v /path/to/data:/data \ - quay.io/cellgeni/adata-cli:latest view /data/your_file.h5ad + quay.io/cellgeni/adata-cli:latest adata view /data/your_file.h5ad ``` ## Documentation diff --git a/pyproject.toml b/pyproject.toml index 1d1ebaf..7767eca 100644 --- a/pyproject.toml +++ b/pyproject.toml @@ -3,7 +3,7 @@ # so the distribution is published as `pyadata-cli`. The import package # and the command are both still `adata`. name = "pyadata-cli" -version = "0.5.0" +version = "0.5.1" description = "Streaming CLI for exploring and editing large AnnData .h5ad and .zarr stores" readme = "README.md" requires-python = ">=3.12" diff --git a/src/adata/storage/__init__.py b/src/adata/storage/__init__.py index 58e7bad..cfea148 100644 --- a/src/adata/storage/__init__.py +++ b/src/adata/storage/__init__.py @@ -440,12 +440,69 @@ def _is_string_src(src: Any) -> bool: return is_string_dtype(getattr(src, "dtype", None)) -def _chunk_step(shape: Sequence[int], chunks: Optional[Sequence[int]]) -> int: - if chunks is not None and len(chunks) > 0 and chunks[0]: - return int(chunks[0]) +#: Byte budget for a single read when streaming a dataset. +#: +#: Using a source's chunk height as the read size means inheriting whatever +#: the writer chose, and a store chunked `(1, n_cols)` is then copied one row +#: per read. On a local disk that is merely wasteful; on a network filesystem +#: (Lustre, NFS) every read is a round-trip costing milliseconds, so a +#: million-row copy spends nearly all of its time waiting. 32 MiB sits +#: comfortably above a typical 1 MiB Lustre stripe while bounding peak memory. +TARGET_READ_BYTES = 32 * 1024 * 1024 + +#: Assumed width of a variable-length string element, for read sizing only. +#: h5py reports `itemsize` 8 for vlen strings because the value is a pointer, +#: which would overestimate the row count by an order of magnitude and break +#: the memory bound. Real cell and gene names sit well under this. +VLEN_ELEMENT_BYTES = 64 + + +def _row_bytes(shape: Sequence[int], dtype: Any) -> int: + """In-memory size of one row along the first axis.""" + width = 1 + for dim in shape[1:]: + width *= max(1, int(dim)) + + itemsize = int(getattr(dtype, "itemsize", 0) or 0) + # 'O' is how h5py spells vlen str, 'T' is numpy StringDType as reported by + # zarr-python 3. Neither itemsize reflects the bytes actually stored. + if itemsize <= 0 or getattr(dtype, "kind", None) in ("O", "T"): + itemsize = VLEN_ELEMENT_BYTES + + return max(1, width * itemsize) + + +def _chunk_step( + shape: Sequence[int], chunks: Optional[Sequence[int]], dtype: Any +) -> int: + """Rows to copy per read, sized for the filesystem rather than the source. + + The step is grown to `TARGET_READ_BYTES`, then rounded down to a whole + number of source chunks: reading part of a chunk still costs decompressing + all of it, so a step that splits one wastes the remainder. + """ if not shape: return 1 - return max(1, min(1024, int(shape[0]))) + + n_rows = int(shape[0]) + if n_rows <= 0: + return 1 + + chunk_rows = 0 + if chunks is not None and len(chunks) > 0 and chunks[0]: + chunk_rows = max(1, int(chunks[0])) + + step = max(1, TARGET_READ_BYTES // _row_bytes(shape, dtype)) + if chunk_rows: + # Never go below a single chunk: a partial read still decompresses the + # whole thing, so a smaller step costs the same I/O for less data. That + # makes TARGET_READ_BYTES a target rather than a cap -- a source whose + # own chunk already exceeds the budget (say 1000 x 1e6 float32, chunked + # whole) reads that chunk regardless. This matches the previous + # behaviour, which used the chunk height verbatim. + step = max(chunk_rows, (step // chunk_rows) * chunk_rows) + + return min(step, n_rows) def copy_dataset(src: Any, dst_group: Any, name: str) -> Any: @@ -475,7 +532,7 @@ def copy_dataset(src: Any, dst_group: Any, name: str) -> Any: ds[()] = src[()] return ds - step = _chunk_step(shape, getattr(src, "chunks", None)) + step = _chunk_step(shape, getattr(src, "chunks", None), src.dtype) for start in range(0, shape[0], step): end = min(start + step, shape[0]) if len(shape) == 1: diff --git a/tests/test_storage.py b/tests/test_storage.py index 926627f..9c57d17 100644 --- a/tests/test_storage.py +++ b/tests/test_storage.py @@ -393,3 +393,77 @@ def test_clamping_handles_one_dimensional_chunks(): assert _clamp_chunks({"chunks": (100,)}, 5)["chunks"] == (5,) assert _clamp_chunks({}, 5) == {} + + +# --------------------------------------------------------------------------- +# read sizing + + +def test_read_step_grows_past_a_pathological_source_chunk(): + """A `(1, n_cols)` source must not be copied one row per read. + + Forwarding the source's chunk height meant inheriting whatever the writer + chose. On a network filesystem each read is a round-trip, so a row-chunked + million-cell store spent the whole copy waiting on latency. + """ + from adata.storage import TARGET_READ_BYTES, _chunk_step, _row_bytes + + shape, chunks, dtype = (1_000_000, 30_000), (1, 30_000), np.dtype("float32") + step = _chunk_step(shape, chunks, dtype) + + assert step > 1, "the source chunk height is no longer used verbatim" + assert step * _row_bytes(shape, dtype) >= TARGET_READ_BYTES // 2 + + +def test_read_step_stays_a_whole_number_of_source_chunks(): + """Reading part of a chunk still costs decompressing all of it.""" + from adata.storage import _chunk_step + + for chunk_rows in (1, 7, 100, 65_536): + step = _chunk_step((10_000_000,), (chunk_rows,), np.dtype("int64")) + assert step % chunk_rows == 0 + + +def test_read_step_never_exceeds_the_dataset(): + from adata.storage import _chunk_step + + assert _chunk_step((10,), (1,), np.dtype("int64")) == 10 + assert _chunk_step((), None, np.dtype("int64")) == 1 + assert _chunk_step((0,), None, np.dtype("int64")) == 1 + + +def test_read_step_bounds_peak_memory(): + """A wide unchunked source used to read 1024 rows however wide they were.""" + from adata.storage import TARGET_READ_BYTES, _chunk_step, _row_bytes + + cases = [ + ((1_000_000, 30_000), None, np.dtype("float32")), + ((1_000_000, 50), (1, 50), np.dtype("float32")), + ((200_000_000,), (65_536,), np.dtype("int64")), + ] + for shape, chunks, dtype in cases: + step = _chunk_step(shape, chunks, dtype) + assert step * _row_bytes(shape, dtype) <= 2 * TARGET_READ_BYTES + + +def test_row_bytes_does_not_trust_a_vlen_itemsize(): + """h5py reports itemsize 8 for vlen str -- that is the pointer, not the text.""" + import h5py + + from adata.storage import VLEN_ELEMENT_BYTES, _row_bytes + + assert _row_bytes((10,), h5py.string_dtype(encoding="utf-8")) == VLEN_ELEMENT_BYTES + assert _row_bytes((10,), np.dtype("O")) == VLEN_ELEMENT_BYTES + assert _row_bytes((10, 4), np.dtype("float32")) == 16 + + +def test_copy_dataset_of_a_row_chunked_source_round_trips(new_store): + """The larger read step must not change what lands on disk.""" + path, opener = new_store() + with opener("a") as root: + from adata.storage import create_dataset + + values = np.arange(400, dtype="float32").reshape(100, 4) + create_dataset(root["uns"], "src", data=values, chunks=(1, 4)) + copy_dataset(root["uns"]["src"], root["uns"], "dst") + assert np.array_equal(root["uns"]["dst"][...], values) diff --git a/uv.lock b/uv.lock index adb6d8d..6c20152 100644 --- a/uv.lock +++ b/uv.lock @@ -486,7 +486,7 @@ wheels = [ [[package]] name = "pyadata-cli" -version = "0.5.0" +version = "0.5.1" source = { editable = "." } dependencies = [ { name = "h5py" },