[Bug] AssertionError: Invalid patching

[jeffquinn-msk]

Description

I have somehow managed to hit an assertion error here:

sopa/sopa/patches/_patches.py

Line 32 in 83789ab

assert (

Command error:
  │ │   patch_width = 300.0                                                    │ │
  │ │         sdata = SpatialData object, with associated Zarr store:          │ │
  │ │                 /data1/tanseyw/projects/nuc2seg/runs/xenium-prime-fresh… │ │
  │ │                 ├── Images                                               │ │
  │ │                 │     └── 'morphology_focus': DataTree[cyx] (5, 101928,  │ │
  │ │                 54013), (5, 50964, 27006), (5, 25482, 13503), (5, 12741, │ │
  │ │                 6751), (5, 6370, 3375)                                   │ │
  │ │                 ├── Points                                               │ │
  │ │                 │     └── 'transcripts': DataFrame with shape:           │ │
  │ │                 (<Delayed>, 14) (3D points)                              │ │
  │ │                 └── Shapes                                               │ │
  │ │                 │     └── 'image_patches': GeoDataFrame shape: (3318, 3) │ │
  │ │                 (2D shapes)                                              │ │
  │ │                 with coordinate systems:                                 │ │
  │ │                 │   ▸ 'global', with elements:                           │ │
  │ │                 │   │   morphology_focus (Images), transcripts (Points), │ │
  │ │                 image_patches (Shapes)                                   │ │
  │ │          self = <repr-error "'OnDiskTranscriptPatches' object has no     │ │
  │ │                 attribute 'bboxes'">                                     │ │
  │ │         tight = True                                                     │ │
  │ │          xmax = np.float32(11475.797)                                    │ │
  │ │          xmin = np.float32(2.671875)                                     │ │
  │ │          ymax = np.float32(21656.938)                                    │ │
  │ │          ymin = np.float32(2.484375)                                     │ │
  │ ╰──────────────────────────────────────────────────────────────────────────╯ │
  │                                                                              │
  │ /opt/bitnami/python/lib/python3.12/site-packages/sopa/patches/_patches.py:33 │
  │ in __init__                                                                  │
  │                                                                              │
  │    30 │   │   if tight:                                                      │
  │    31 │   │   │   self.patch_width = self.tight_width()                      │
  │    32 │   │   │   assert (                                                   │
  │ ❱  33 │   │   │   │   self._count == self.count()                            │
  │    34 │   │   │   ), f"Invalid patching with {self.delta=}, {self.patch_widt │
  │    35 │                                                                      │
  │    36 │   def count(self):                                                   │
  │                                                                              │
  │ ╭───────────────────────────────── locals ─────────────────────────────────╮ │
  │ │    int_coords = False                                                    │ │
  │ │ patch_overlap = 30.0                                                     │ │
  │ │   patch_width = 300.0                                                    │ │
  │ │          self = <sopa.patches._patches.Patches1D object at               │ │
  │ │                 0x7f6c5c006f00>                                          │ │
  │ │         tight = True                                                     │ │
  │ │          xmax = np.float32(11475.797)                                    │ │
  │ │          xmin = np.float32(2.671875)                                     │ │
  │ ╰──────────────────────────────────────────────────────────────────────────╯ │
  ╰──────────────────────────────────────────────────────────────────────────────╯
  AssertionError: Invalid patching with self.delta=np.float32(11473.125),
  self.patch_width=np.float32(296.11917) and self.patch_overlap=30.0

Reproducing the issue

Apptainer> python
Python 3.12.8 (main, Jan 17 2025, 18:22:41) [GCC 12.2.0] on linux
Type "help", "copyright", "credits" or "license" for more information.
>>> import numpy as np
>>> from sopa.patches._patches import Patches1D
>>> xmax = np.float32(11475.797)
>>> xmin = np.float32(2.671875)
>>> tight =  True
>>> patch_width = 300.0
>>> patch_overlap = 30.0
>>> int_coords = False
>>> Patches1D(xmin=xmin, xmax=xmax, patch_width=patch_width, patch_overlap=patch_overlap, tight=tight, int_coords=int_coords)
Traceback (most recent call last):
  File "<stdin>", line 1, in <module>
  File "/opt/bitnami/python/lib/python3.12/site-packages/sopa/patches/_patches.py", line 33, in __init__
    self._count == self.count()
AssertionError: Invalid patching with self.delta=np.float32(11473.125), self.patch_width=np.float32(296.11917) and self.patch_overlap=30.0

Note: I was not able to reproduce the error on Python 3.11.11, so that gives some hints

Expected behavior

I've read over it a bit but I'm still not clear what this assertion error is checking for. I would expect it not to fail given that my patching parameters seem reasonable, the overlap is smaller than the patch width, which is smaller than xmax-xmin.

Weirder still, a small perturbation in the params results in no error:

>>> Patches1D(xmin=xmin, xmax=xmax, patch_width=patch_width, patch_overlap=patch_overlap-1, tight=tight, int_coords=int_coords)
<sopa.patches._patches.Patches1D object at 0x7f8dee2a2bd0>

System

Apptainer> python --version
Python 3.12.8
Apptainer> cat /etc/os-release
PRETTY_NAME="Debian GNU/Linux 12 (bookworm)"
NAME="Debian GNU/Linux"
VERSION_ID="12"
VERSION="12 (bookworm)"
VERSION_CODENAME=bookworm
ID=debian
HOME_URL="https://www.debian.org/"
SUPPORT_URL="https://www.debian.org/support"
BUG_REPORT_URL="https://bugs.debian.org/"
Apptainer> pip freeze
aiobotocore==2.19.0
aiohappyeyeballs==2.4.6
aiohttp==3.11.12
aioitertools==0.12.0
aiosignal==1.3.2
anndata==0.11.3
annotated-types==0.7.0
array_api_compat==1.10.0
asciitree==0.3.3
attrs==25.1.0
botocore==1.36.3
cellpose==3.1.0
certifi==2025.1.31
charset-normalizer==3.4.1
click==8.1.8
cloudpickle==3.1.1
colorcet==3.1.0
contourpy==1.3.1
cycler==0.12.1
dask==2024.11.2
dask-expr==1.1.19
dask-image==2024.5.3
datashader==0.17.0
Deprecated==1.2.18
distlib==0.3.9
distributed==2024.11.2
fasteners==0.19
fastremap==1.15.0
filelock==3.16.1
flake8==7.1.1
FlowIO==1.3.0
fonttools==4.56.0
frozenlist==1.5.0
fsspec==2024.12.0
geopandas==1.0.1
h5py==3.12.1
idna==3.10
imagecodecs==2024.12.30
imageio==2.37.0
Jinja2==3.1.5
jmespath==1.0.1
joblib==1.4.2
kiwisolver==1.4.8
lazy_loader==0.4
legacy-api-wrap==1.4.1
llvmlite==0.44.0
locket==1.0.0
loompy==3.0.8
markdown-it-py==3.0.0
MarkupSafe==3.0.2
matplotlib==3.10.0
mccabe==0.7.0
mdurl==0.1.2
mpmath==1.3.0
msgpack==1.1.0
multidict==6.1.0
multipledispatch==1.0.0
multiscale_spatial_image==2.0.2
natsort==8.4.0
networkx==3.4.2
numba==0.61.0
numcodecs==0.15.1
numpy==2.0.2
numpy-groupies==0.11.2
nvidia-cublas-cu12==12.4.5.8
nvidia-cuda-cupti-cu12==12.4.127
nvidia-cuda-nvrtc-cu12==12.4.127
nvidia-cuda-runtime-cu12==12.4.127
nvidia-cudnn-cu12==9.1.0.70
nvidia-cufft-cu12==11.2.1.3
nvidia-curand-cu12==10.3.5.147
nvidia-cusolver-cu12==11.6.1.9
nvidia-cusparse-cu12==12.3.1.170
nvidia-cusparselt-cu12==0.6.2
nvidia-nccl-cu12==2.21.5
nvidia-nvjitlink-cu12==12.4.127
nvidia-nvtx-cu12==12.4.127
ome-types==0.5.3
ome-zarr==0.10.3
opencv-python==4.11.0.86
opencv-python-headless==4.11.0.86
packaging==24.2
pandas==2.2.3
param==2.2.0
partd==1.4.2
patsy==1.0.1
pillow==11.1.0
PIMS==0.7
platformdirs==4.3.6
pooch==1.8.2
propcache==0.2.1
psutil==6.1.1
pyarrow==19.0.0
pycodestyle==2.12.1
pyct==0.5.0
pydantic==2.10.6
pydantic-compat==0.1.2
pydantic_core==2.27.2
pyflakes==3.2.0
Pygments==2.19.1
pynndescent==0.5.13
pyogrio==0.10.0
pyparsing==3.2.1
pyproj==3.7.0
python-dateutil==2.9.0.post0
pytz==2025.1
PyYAML==6.0.2
readfcs==2.0.1
requests==2.32.3
rich==13.9.4
roifile==2024.9.15
s3fs==2024.12.0
scanpy==1.10.4
scikit-image==0.25.1
scikit-learn==1.6.1
scipy==1.15.1
seaborn==0.13.2
session_info==1.0.0
setuptools @ file:///bitnami/blacksmith-sandox/setuptools-70.3.0
shapely==2.0.7
shellingham==1.5.4
six==1.17.0
slicerator==1.1.0
sortedcontainers==2.4.0
spatial_image==1.1.0
spatialdata==0.3.0
spatialdata-io==0.1.7
statsmodels==0.14.4
stdlib-list==0.11.0
sympy==1.13.1
tblib==3.0.0
threadpoolctl==3.5.0
tifffile==2025.1.10
toml==0.10.2
toolz==1.0.0
torch==2.6.0
torchaudio==2.6.0
torchvision==0.21.0
tornado==6.4.2
tqdm==4.67.1
triton==3.2.0
typer==0.15.1
typing_extensions==4.12.2
tzdata==2025.1
umap-learn==0.5.7
urllib3==2.3.0
virtualenv @ file:///bitnami/blacksmith-sandox/virtualenv-20.29.1
wrapt==1.17.2
xarray==2024.11.0
xarray-dataclasses==1.9.1
xarray-schema==0.0.3
xarray-spatial==0.4.0
xsdata==24.3.1
yarl==1.18.3
zarr==2.18.4
zict==3.0.0

[quentinblampey]

Hi @jeffquinn-msk, I was able to reproduce the issue. For some context: The Patches1D class tries to use the smallest patch_width that keeps the same number of patches along the axis.

This is because, if you have a width of 1000 and a patch-width of 950, you don't want one patch of size 950 and one patch of size 50. Instead, we prefer having two patches of size 500 (here, for this example, I assume there is no overlap). Then, the purpose of the assert is to make sure that, after defining the "tight patch width," we still have the same number of patches along the axis, which is, by definition, what we want to achieve.

Note that we update the patch_width as below, where EPS=1e-5 is used exactly to avoid the rounding issues that you experience.

patch_width = (self.delta + (self._count - 1) * self.patch_overlap) / self._count + EPS

In your case, it seems that EPS is a little bit too low. Using EPS=1e-4 solves the issue. It's weird that this issue depends on the version of Python...And what is also weird is that if you use normal float values instead of np.float32, it also solves the issue:

xmax = 11475.797 # instead of np.float32(11475.797)
xmin = 2.671875 # instead of np.float32(2.671875)

Anyway, I can update EPS, but I think this is really an edge case, probably due to the combination of multiple factors: (i) having a very large delta compared to a very small patch_width, (ii) using np.float32 for xmax / xmin, and (iii) maybe some bad luck for the parameter choice.

[jeffquinn-msk]

And what is also weird is that if you use normal float values instead of np.float32, it also solves the issue:

Yes I noticed that as well.

(i) having a very large delta compared to a very small patch_width

These are reasonable parameters in my experience for running Baysor on Xenium 5k datasets

(iii) maybe some bad luck for the parameter choice.

My luck must be very bad indeed, as I have repeatedly hit this problem on different datasets, each time just permuting the patch size a bit.

[quentinblampey]

These are reasonable parameters in my experience for running Baysor on Xenium 5k datasets

May I ask how much RAM it uses for each patch?

My luck must be very bad indeed, as I have repeatedly hit this problem on different datasets, each time just permuting the patch size a bit.

Okay I see. So I pushed a change in EPS to prevent this.
Let me know if the issue happens again - in that case I would try to cast xmax and xmin to float values inside the Patch1D class.

[jeffquinn-msk]

May I ask how much RAM it uses for each patch?

These were the memory usage characteristics for one of the recent runs with this patch size. It is spiky, given that it depends also on the number of transcripts in each patch and not just the patch size, hence why I had to set the patch size smaller to avoid OOMs. Runtime also is scaling non-linearly with patch size which is another reason why I resorted to the small patches.

Okay I see. So I pushed a change in EPS to prevent this.

Ok sounds good! The root of this seems to be a floating point precision problem, but I don't understand quite why its happening, will need to sit down and stare at the code a bit more.

It makes sense that float32 would make it worse in that case, given that a builtin python float has 64 bits. Also there were some floating point arithmetic changes in 3.12 vs 3.11, such as the sum() implementation for floats, that could explain the version dependent behavior:

https://docs.python.org/3/whatsnew/3.12.html

[quentinblampey]

Thanks for the detailed memory summary!

The root of this seems to be a floating point precision problem, but I don't understand quite why its happening

Yes, it's definitely a floating point precision issue. I probably won't spend a lot of time digging into this, but if you do please let me know I'm curious about it!

I'm closing the issue as it will be solved in sopa==2.0.2, which is expected to be released either tomorrow or next week.

[quentinblampey]

The new version is released!