Visium hd
This adds support for annotating the resulting tables with corresponding labels layers (the bins).
Annotating the table by a labels layer has the advantage that the bins, and their expression levels can be viewed in napari-spatialdata, which is not possible if you annotate via the shapes layers ( it is possible when setting bins_as_squares to False, but this will give you circles, not squares, which is not ideal ).
Viewing and loading expression levels for 016um and 008um bins is quite fast, viewing the 002um bins is also still possible, althoug loading expression levels associated with the bins takes some seconds for each gene, but once loaded is quite fast and convenient to view.
See screenshot for 008um:
It is not as snappy as https://github.com/scverse/spatialdata/blob/8879aff17b8169c6bb6ff3537e4ddcee0204f168/src/spatialdata/_core/operations/rasterize_bins.py#L28, but maybe a little bit more convenient in use.
If this PR would be merged, then this function spatialdata.rasterize_bins would also need to be updated. I would remove the bins parameter, and estimate the transformation only using the tables, i.e. using table.obsm[ 'spatial' ] and the col_key and row_key in table.obs; or get the transformation from the labels layer.
Codecov Report
Attention: Patch coverage is 23.07692% with 10 lines in your changes missing coverage. Please review.
Project coverage is 46.44%. Comparing base (
c4f6cf7) to head (261dbbb). Report is 135 commits behind head on main.
| Files with missing lines | Patch % | Lines |
|---|---|---|
| src/spatialdata_io/readers/visium_hd.py | 23.07% | 10 Missing :warning: |
Additional details and impacted files
@@ Coverage Diff @@
## main #211 +/- ##
==========================================
+ Coverage 43.68% 46.44% +2.76%
==========================================
Files 23 23
Lines 2353 2394 +41
==========================================
+ Hits 1028 1112 +84
+ Misses 1325 1282 -43
| Files with missing lines | Coverage Δ | |
|---|---|---|
| src/spatialdata_io/readers/visium_hd.py | 18.28% <23.07%> (-0.01%) |
:arrow_down: |
🚀 New features to boost your workflow:
- ❄ Test Analytics: Detect flaky tests, report on failures, and find test suite problems.
Thanks for the PR, it looks great, very curios to experiment with it! I have some tasks on my todo list but I hope to get to this soon.
Hi, great PR, it looks very interesting! I´m trying to use it for my project, but ran into an error. I´m relatively new to using PRs on github, so maybe I already did a mistake during setup. I created a new conda environment and installed your branch via "pip install git+https://github.com/ArneDefauw/spatialdata-io.git@visium_hd", added an ipykernel for jupyter and tried to use the visium_hd reader function in a notebook. I specified my paths, but get a FileNotFoundError for following code:
from spatialdata_io import visium_hd sdata = visium_hd( path_read, load_all_images=True, annotate_table_by_labels=True, fullres_image_file=fullres_image_file, dataset_id=dataset_id )
FileNotFoundError: [Errno 2] Unable to synchronously open file (unable to open file: name = '/files/03-SpaceRanger-Results/P31904_402/outs/P31904_402_feature_slice.h5', errno = 2, error message = 'No such file or directory', flags = 0, o_flags = 0)
"P31904_402_feature_slice.h5" indeed does not exist, it should simply be "feature_slice.h5" in "outs" of the "P31904_402" folder. I don´t understand why the error is thrown though, since that part of the code was not changed to my understanding?
Here is the full traceback: FileNotFoundError Traceback (most recent call last) Cell In[4], line 2 1 from spatialdata_io import visium_hd ----> 2 sdata = visium_hd( 3 path_read, 4 load_all_images=True, 5 annotate_table_by_labels=True, 6 fullres_image_file=fullres_image_file, 7 dataset_id=dataset_id 8 )
File ~/anaconda3/envs/sd_pr_test/lib/python3.12/site-packages/spatialdata_io/readers/visium_hd.py:115, in visium_hd(path, dataset_id, filtered_counts_file, bin_size, bins_as_squares, annotate_table_by_labels, fullres_image_file, load_all_images, imread_kwargs, image_models_kwargs, anndata_kwargs) 104 def load_image(path: Path, suffix: str, scale_factors: list[int] | None = None) -> None: 105 _load_image( 106 path=path, 107 images=images, (...) 112 scale_factors=scale_factors, 113 ) --> 115 metadata, hd_layout = _parse_metadata(path, filename_prefix) 116 transform_matrices = _get_transform_matrices(metadata, hd_layout) 117 file_format = hd_layout[VisiumHDKeys.FILE_FORMAT]
File ~/anaconda3/envs/sd_pr_test/lib/python3.12/site-packages/spatialdata_io/readers/visium_hd.py:471, in _parse_metadata(path, filename_prefix) 470 def _parse_metadata(path: Path, filename_prefix: str) -> tuple[dict[str, Any], dict[str, Any]]: --> 471 with h5py.File(path / f"{filename_prefix}{VisiumHDKeys.FEATURE_SLICE_FILE.value}", "r") as f5: 472 metadata = json.loads(dict(f5.attrs)[VisiumHDKeys.METADATA_JSON]) 473 hd_layout = json.loads(metadata[VisiumHDKeys.HD_LAYOUT_JSON])
File ~/anaconda3/envs/sd_pr_test/lib/python3.12/site-packages/h5py/_hl/files.py:562, in File.init(self, name, mode, driver, libver, userblock_size, swmr, rdcc_nslots, rdcc_nbytes, rdcc_w0, track_order, fs_strategy, fs_persist, fs_threshold, fs_page_size, page_buf_size, min_meta_keep, min_raw_keep, locking, alignment_threshold, alignment_interval, meta_block_size, **kwds) 553 fapl = make_fapl(driver, libver, rdcc_nslots, rdcc_nbytes, rdcc_w0, 554 locking, page_buf_size, min_meta_keep, min_raw_keep, 555 alignment_threshold=alignment_threshold, 556 alignment_interval=alignment_interval, 557 meta_block_size=meta_block_size, 558 **kwds) 559 fcpl = make_fcpl(track_order=track_order, fs_strategy=fs_strategy, 560 fs_persist=fs_persist, fs_threshold=fs_threshold, 561 fs_page_size=fs_page_size) --> 562 fid = make_fid(name, mode, userblock_size, fapl, fcpl, swmr=swmr) 564 if isinstance(libver, tuple): 565 self._libver = libver
File ~/anaconda3/envs/sd_pr_test/lib/python3.12/site-packages/h5py/_hl/files.py:235, in make_fid(name, mode, userblock_size, fapl, fcpl, swmr) 233 if swmr and swmr_support: 234 flags |= h5f.ACC_SWMR_READ --> 235 fid = h5f.open(name, flags, fapl=fapl) 236 elif mode == 'r+': 237 fid = h5f.open(name, h5f.ACC_RDWR, fapl=fapl)
File h5py/_objects.pyx:54, in h5py._objects.with_phil.wrapper()
File h5py/_objects.pyx:55, in h5py._objects.with_phil.wrapper()
File h5py/h5f.pyx:102, in h5py.h5f.open()
Hi, great PR, it looks very interesting! I´m trying to use it for my project, but ran into an error. I´m relatively new to using PRs on github, so maybe I already did a mistake during setup. I created a new conda environment and installed your branch via "pip install git+https://github.com/ArneDefauw/spatialdata-io.git@visium_hd", added an ipykernel for jupyter and tried to use the visium_hd reader function in a notebook. I specified my paths, but get a FileNotFoundError for following code:
from spatialdata_io import visium_hd sdata = visium_hd( path_read, load_all_images=True, annotate_table_by_labels=True, fullres_image_file=fullres_image_file, dataset_id=dataset_id )
FileNotFoundError: [Errno 2] Unable to synchronously open file (unable to open file: name = '/files/03-SpaceRanger-Results/P31904_402/outs/P31904_402_feature_slice.h5', errno = 2, error message = 'No such file or directory', flags = 0, o_flags = 0)
"P31904_402_feature_slice.h5" indeed does not exist, it should simply be "feature_slice.h5" in "outs" of the "P31904_402" folder. I don´t understand why the error is thrown though, since that part of the code was not changed to my understanding?
Here is the full traceback: FileNotFoundError Traceback (most recent call last) Cell In[4], line 2 1 from spatialdata_io import visium_hd ----> 2 sdata = visium_hd( 3 path_read, 4 load_all_images=True, 5 annotate_table_by_labels=True, 6 fullres_image_file=fullres_image_file, 7 dataset_id=dataset_id 8 )
File ~/anaconda3/envs/sd_pr_test/lib/python3.12/site-packages/spatialdata_io/readers/visium_hd.py:115, in visium_hd(path, dataset_id, filtered_counts_file, bin_size, bins_as_squares, annotate_table_by_labels, fullres_image_file, load_all_images, imread_kwargs, image_models_kwargs, anndata_kwargs) 104 def load_image(path: Path, suffix: str, scale_factors: list[int] | None = None) -> None: 105 _load_image( 106 path=path, 107 images=images, (...) 112 scale_factors=scale_factors, 113 ) --> 115 metadata, hd_layout = _parse_metadata(path, filename_prefix) 116 transform_matrices = _get_transform_matrices(metadata, hd_layout) 117 file_format = hd_layout[VisiumHDKeys.FILE_FORMAT]
File ~/anaconda3/envs/sd_pr_test/lib/python3.12/site-packages/spatialdata_io/readers/visium_hd.py:471, in _parse_metadata(path, filename_prefix) 470 def _parse_metadata(path: Path, filename_prefix: str) -> tuple[dict[str, Any], dict[str, Any]]: --> 471 with h5py.File(path / f"{filename_prefix}{VisiumHDKeys.FEATURE_SLICE_FILE.value}", "r") as f5: 472 metadata = json.loads(dict(f5.attrs)[VisiumHDKeys.METADATA_JSON]) 473 hd_layout = json.loads(metadata[VisiumHDKeys.HD_LAYOUT_JSON])
File ~/anaconda3/envs/sd_pr_test/lib/python3.12/site-packages/h5py/_hl/files.py:562, in File.init(self, name, mode, driver, libver, userblock_size, swmr, rdcc_nslots, rdcc_nbytes, rdcc_w0, track_order, fs_strategy, fs_persist, fs_threshold, fs_page_size, page_buf_size, min_meta_keep, min_raw_keep, locking, alignment_threshold, alignment_interval, meta_block_size, **kwds) 553 fapl = make_fapl(driver, libver, rdcc_nslots, rdcc_nbytes, rdcc_w0, 554 locking, page_buf_size, min_meta_keep, min_raw_keep, 555 alignment_threshold=alignment_threshold, 556 alignment_interval=alignment_interval, 557 meta_block_size=meta_block_size, 558 **kwds) 559 fcpl = make_fcpl(track_order=track_order, fs_strategy=fs_strategy, 560 fs_persist=fs_persist, fs_threshold=fs_threshold, 561 fs_page_size=fs_page_size) --> 562 fid = make_fid(name, mode, userblock_size, fapl, fcpl, swmr=swmr) 564 if isinstance(libver, tuple): 565 self._libver = libver
File ~/anaconda3/envs/sd_pr_test/lib/python3.12/site-packages/h5py/_hl/files.py:235, in make_fid(name, mode, userblock_size, fapl, fcpl, swmr) 233 if swmr and swmr_support: 234 flags |= h5f.ACC_SWMR_READ --> 235 fid = h5f.open(name, flags, fapl=fapl) 236 elif mode == 'r+': 237 fid = h5f.open(name, h5f.ACC_RDWR, fapl=fapl)
File h5py/_objects.pyx:54, in h5py._objects.with_phil.wrapper()
File h5py/_objects.pyx:55, in h5py._objects.with_phil.wrapper()
File h5py/h5f.pyx:102, in h5py.h5f.open()
Hi @t-a-m-i ,
Does this work using spatialdata-io (last release), installed via pip install spatialdata-io? And then do:
from spatialdata_io import visium_hd
sdata = visium_hd(
path_read,
load_all_images=True,
fullres_image_file=fullres_image_file,
dataset_id=dataset_id
)
I don't think the error you encounter is related to this PR.
yes, I was also confused about the error since it didn´t seem related to the PR. I didn´t realise though that there was a new release. I get the same error with the latest version of spatialdata-io (v0.1.5), with v0.1.4 everything works fine, I´ll do a bug report then, thanks!
Thanks @t-a-m-i for opening https://github.com/scverse/spatialdata-io/issues/216, I'll follow up there.
Notebook relevant for the PR review: https://scverse.zulipchat.com/#narrow/channel/315824-spatial/topic/SpatialData/near/479568220
Thanks again @ArneDefauw and @t-a-m-i. Today I studied your code and reviewed and tested this PR. Here are my comments.
-
The PR works great and adds an important use case, because unlike the on-demand
rasterize_bins(), here we can save the data to disk (as observed the performance is slower than the on-demandrasterize_bins(), but faster than when plotting shapes). -
The example from @t-a-m-i shows that sometimes one may want to convert the shapes from an existing SpatialData
.zarrVisium HD object to labels. The approach that she took was to extend therasterize_bins()functions to replicate the code introduced in this PR. For this reason I believe that the code in this PR would better suited directly inspatialdatarather than inspatialdata-io. Some further arguments for this are that in this way this function would be accessible also for other technologies, such as Stereo-seq and Open-ST, and the code would be cleaner since right now some parts duplicate code already present inrasterize_bins().I therefore kindly ask you @ArneDefauw if you could open a new PR in
spatialdataand move the current code there; in the current PR you would then simply call the new implementation ofrasterize_bins(). To distinguish between the on-demandrasterize_bins()(which returns images) and therasterize_bins()returning labels, I would add a new argumentreturn_regions_as_labels: bool = False(as we do forrasterize()). Also, inspatialdatait would be important to include tests for this new functionality. -
Last comment, there is only one logic in the code that I'd change: I would not shift the instance_id by 1 because this could break compatibility (or create confusion) with previously written datasets; rather I would add a new column to the table to be used only when mapping the table to the newly created labels object.
Hi @LucaMarconato , thanks for having a look! It makes sense to move the code for this PR to spatialdata, and I probably have some time next week to do it. I will then also update this PR.
I think this PR is ready for review again, together with related PR https://github.com/scverse/spatialdata/pull/811
Note, we probably need to update this notebook https://spatialdata.scverse.org/en/stable/tutorials/notebooks/notebooks/examples/technology_visium_hd.html after the PR's are merged. I would be happy to help.
Thanks @ArneDefauw for updating this PR. I have reviewed it and it's ready to merge. Please notice that I have reverted the change on the instance_key column 8c39936 (#211). Otherwise we could break backward compatibility with scripts/notebooks built with previous versions of spatialdata-io.