Chunkwise image loader

[lucas-diedrich]

Description

This PR addresses the challenge that the currently implemented and planned image loaders require loading imaging data entirely into memory, typically as NumPy arrays. Given the large size of microscopy datasets, this is not always feasible.

To mitigate this issue, and as discussed with @LucaMarconato, this PR aims to introduce a generalizable approach for reading large microscopy files in chunks, enabling efficient handling of data that does not fit into memory.

Some related discussions.

• Image Functions spatialdata#704 (comment)
• Efficient whole slide imaging IO spatialdata#856 (comment)
• DVP image readers lucas-diedrich/spatialdata-io#1 (comment)

Strategy

In this PR, we focus on .tiff images, as implemented in the _tiff_to_chunks function.

1. Get a lazy representation of the image via a suitable reader function (here: tifffile.memmap)
2. Pre-define chunks that fit into memory, based on the dimensions of the image (_compute_chunks)
3. Load small chunks via a custom reader function and pass the chunks to dask.array which is memory-mapped and avoids memory overflow (_read_chunks)
4. Reassembling the chunks into a dask.array (via dask.array.block)
5. Parse to Image2DModel.

The strategy is implemented in

• src/spatialdata_io/readers/generic.py and
• src/spatialdata_io/readers/_utils/_image.py

Future extensions

The strategy can be implemented for any image type, as long as it is possible to implement

1. a lazy image-data loader
2. define a custom reader function

We have implemented similar readers for openslide-compatible whole slide images and the Carl-Zeiss microscopy format.

[codecov-commenter]

Codecov Report

❌ Patch coverage is 95.69892% with 4 lines in your changes missing coverage. Please review.
✅ Project coverage is 62.91%. Comparing base (2ebfab0) to head (4b4a00e).
⚠️ Report is 16 commits behind head on main.

Files with missing lines	Patch %	Lines
src/spatialdata_io/readers/generic.py	91.11%	4 Missing ⚠️

Additional details and impacted files

@@            Coverage Diff             @@
##             main     #279      +/-   ##
==========================================
+ Coverage   55.16%   62.91%   +7.74%     
==========================================
  Files          26       27       +1     
  Lines        2844     3117     +273     
==========================================
+ Hits         1569     1961     +392     
+ Misses       1275     1156     -119     

Files with missing lines	Coverage Δ
src/spatialdata_io/readers/_utils/_image.py	100.00% <100.00%> (ø)
src/spatialdata_io/readers/generic.py	88.88% <91.11%> (+3.17%)	⬆️

... and 4 files with indirect coverage changes

🚀 New features to boost your workflow:

• ❄️ Test Analytics: Detect flaky tests, report on failures, and find test suite problems.

[melonora]

Thanks for your contribution! I have 2 minor suggestions. I also saw that you use the width by height convention. Personally, I don't have a strong opinion here, though we could also stick to array api conventions. @LucaMarconato WDYT? Pre-approving for now.

[melonora]

sorry had to change due to rethinking memmap. This does not always work, for example when dealing with compressed tiffs as far as I am aware.

[lucas-diedrich]

Thanks! Addressed your comments

• Naming, simplification 46ed3e5, 99731fe
• Read compressed images with dask_image.imread instead of tifffile 9e057de

[lucas-diedrich]

Hi @LucaMarconato - Just saw that you pushed some updates. Could you comment on the current implementation/is this PR still interesting for you?

[LucaMarconato]

@lucas-diedrich yes still interested in it, I have time now to go through it. I added my "code review" in terms of "TODO" items for myself. I'll check the code and push some modifications if I see that some minor fixes are needed; I'll comment eventual larger changes, but the code looks good. I'll also do a benchmark with asv for performance.

[LucaMarconato]

There was a problem with the dimension of the returned data. Image2DModel.parse() always returns a cyx, but the variable reference in test_read_tiff() was showing a non-cyx image. Now it's fixed.

[LucaMarconato]

@lucas-diedrich I'm done with edits, please double check the changes if you have time. The main comment, added in the notes in _read_chunks(), is the following:

    Notes
    -------
    As shown in `_compute_chunks()`, `coords` are in the form `(x, y,
    width, height)`. In that function, the inner list (dim = -1) iterates over `y`
    values, and the outer list (dim = -2) iterates over `x` values. In `_read_chunks(
    )`, we use the more common `(y, x)` ordering: the inner list (dim = -1) iterates
    over `x` values, and the outer list (dim = -2) iterates over `y` values.

    This mismatch can be confusing. A straightforward fix is to standardize `coords`
    to `(y, x, height, width)` instead of `(x, y, width, height)`.

In summary, redefining coords as y, x, height, width would make the code more readable, but it works also as it is, so no strong opinion.

I will now do the benchmark with asv.

[lucas-diedrich]

Hi @LucaMarconato, thanks for implementing the changes and apologies for the confusing coords convention!

The tests are currently failing when I pass an asymmetric chunk size (e.g (29, 71), ((30, 50) in test_read_tiff), I presume as there is a subtle error with the dimension ordering. I'll look into it.