Merge pull request #45 from gustaveroussy/dev

hotfix v1.0.8

.github/workflows/ci.yml

@@ -58,14 +58,17 @@ jobs:
       #----------------------------------------------
       - name: Install dependencies
         if: steps.cached-poetry-dependencies.outputs.cache-hit != 'true'
-        run: poetry install --extras "dev cellpose" --no-interaction
+        run: poetry install --extras "dev cellpose snakemake" --no-interaction
 
       #----------------------------------------------
       #              perform tasks
       #----------------------------------------------
       - name: Tests
         run: poetry run pytest
 
+      - name: Snakemake
+        run: cd workflow && poetry run snakemake --config sdata_path=tuto.zarr --configfile=config/toy/uniform_cellpose.yaml -c1
+
       - name: Deploy doc
         if: contains(github.ref, 'tags')
         run: poetry run mkdocs gh-deploy --force

CHANGELOG.md

@@ -1,3 +1,7 @@
+## [1.0.8] - 2024-04-02
+
+Hotfix: resolve issues related to `spatialdata>=1.0.0`
+
 ## [1.0.7] - 2024-03-29
 
 ### Changed

pyproject.toml

@@ -1,6 +1,6 @@
 [tool.poetry]
 name = "sopa"
-version = "1.0.7"
+version = "1.0.8"
 description = "Spatial-omics pipeline and analysis"
 documentation = "https://gustaveroussy.github.io/sopa"
 homepage = "https://gustaveroussy.github.io/sopa"
@@ -23,8 +23,8 @@ sopa = "sopa.main:app"
 
 [tool.poetry.dependencies]
 python = ">=3.9,<3.12"
-spatialdata = ">=0.1.1"
-spatialdata-io = ">=0.1.1"
+spatialdata = ">=0.1.2"
+spatialdata-io = ">=0.1.2"
 scanpy = ">=1.9.5,!=1.9.7"
 botocore = "1.34.19"
 typer = ">=0.9.0"

sopa/_sdata.py

@@ -1,14 +1,18 @@
 from __future__ import annotations
 
 import logging
+from pathlib import Path
 from typing import Iterator
 
 import geopandas as gpd
 import pandas as pd
 import xarray as xr
+import zarr
 from multiscale_spatial_image import MultiscaleSpatialImage
+from ome_zarr.io import parse_url
 from spatial_image import SpatialImage
 from spatialdata import SpatialData
+from spatialdata._io import write_image, write_shapes, write_table
 from spatialdata.models import SpatialElement
 from spatialdata.transformations import Identity, get_transformation, set_transformation
 
@@ -176,3 +180,44 @@ def get_spatial_image(
     if return_key:
         return key, image
     return image
+
+
+def save_shapes(
+    sdata: SpatialData,
+    name: str,
+    overwrite: bool = False,
+) -> None:
+    elem_group = sdata._init_add_element(name=name, element_type="shapes", overwrite=overwrite)
+    write_shapes(
+        shapes=sdata.shapes[name],
+        group=elem_group,
+        name=name,
+    )
+
+
+def save_image(
+    sdata: SpatialData,
+    name: str,
+    overwrite: bool = False,
+) -> None:
+    elem_group = sdata._init_add_element(name=name, element_type="images", overwrite=overwrite)
+    write_image(
+        image=sdata.images[name],
+        group=elem_group,
+        name=name,
+    )
+    from spatialdata._io.io_raster import _read_multiscale
+
+    # reload the image from the Zarr storage so that now the element is lazy loaded, and most importantly,
+    # from the correct storage
+    assert elem_group.path == "images"
+    path = Path(elem_group.store.path) / "images" / name
+    image = _read_multiscale(path, raster_type="image")
+    sdata._add_image_in_memory(name=name, image=image, overwrite=True)
+
+
+def save_table(sdata: SpatialData):
+    store = parse_url(sdata.path, mode="r+").store
+    root = zarr.group(store=store)
+    elem_group = root.require_group(name="table")
+    write_table(table=sdata.table, group=elem_group, name="table")

sopa/embedding/patches.py

@@ -17,10 +17,10 @@
 from spatialdata.models import Image2DModel
 from spatialdata.transformations import Scale
 
-import sopa.embedding.models as models
-from sopa._constants import SopaKeys
-from sopa._sdata import get_intrinsic_cs, get_key
-from sopa.segmentation import Patches2D
+from .._constants import SopaKeys
+from .._sdata import get_intrinsic_cs, get_key, save_image
+from ..segmentation import Patches2D
+from . import models
 
 log = logging.getLogger(__name__)
 
@@ -217,7 +217,8 @@ def embed_wsi_patches(
     embedding_image.coords["x"] = embedder.patch_width * embedding_image.coords["x"]
 
     embedding_key = f"sopa_{model_name}"
-    sdata.add_image(embedding_key, embedding_image)
+    sdata.images[embedding_key] = embedding_image
+    save_image(sdata, embedding_key)
 
     log.info(f"WSI embeddings saved as an image in sdata['{embedding_key}']")
 

sopa/io/explorer/images.py

@@ -13,7 +13,7 @@
 from spatialdata.transformations import Affine, set_transformation
 from tqdm import tqdm
 
-from ..._sdata import get_intrinsic_cs, get_spatial_image
+from ..._sdata import get_intrinsic_cs, get_spatial_image, save_image
 from ...utils.image import resize_numpy, scale_dtype
 from ._constants import ExplorerConstants, FileNames, image_metadata
 from .utils import explorer_file_path
@@ -224,4 +224,5 @@ def align(
     set_transformation(image, {pixel_cs: to_pixel}, set_all=True)
 
     log.info(f"Adding image {image_name}:\n{image}")
-    sdata.add_image(image_name, image, overwrite=overwrite)
+    sdata.images[image_name] = image
+    save_image(sdata, image_name, overwrite=overwrite)

sopa/segmentation/aggregate.py

@@ -25,6 +25,8 @@
     get_element,
     get_item,
     get_spatial_image,
+    save_shapes,
+    save_table,
     to_intrinsic,
 )
 from ..io.explorer.utils import str_cell_id
@@ -73,7 +75,8 @@ def standardize_table(self):
         self.table.obs_names = list(map(str_cell_id, range(self.table.n_obs)))
 
         self.geo_df.index = list(self.table.obs_names)
-        self.sdata.add_shapes(self.shapes_key, self.geo_df, overwrite=True)
+        self.sdata.shapes[self.shapes_key] = self.geo_df
+        save_shapes(self.sdata, self.shapes_key, overwrite=True)
 
         self.table.obsm["spatial"] = np.array(
             [[centroid.x, centroid.y] for centroid in self.geo_df.centroid]
@@ -100,7 +103,8 @@ def standardize_table(self):
 
         if self.sdata.table is not None and self.overwrite:
             del self.sdata.table
-        self.sdata.table = self.table
+        self.sdata.tables["table"] = self.table
+        save_table(self.sdata)
 
     def filter_cells(self, where_filter: np.ndarray):
         log.info(f"Filtering {where_filter.sum()} cells")

sopa/segmentation/baysor/resolve.py

@@ -16,7 +16,7 @@
 from tqdm import tqdm
 
 from ..._constants import SopaKeys
-from ..._sdata import get_element, get_key
+from ..._sdata import get_element, get_key, save_shapes
 from .. import aggregate, shapes
 
 log = logging.getLogger(__name__)
@@ -171,7 +171,8 @@ def resolve(
         instance_key=SopaKeys.INSTANCE_KEY,
     )
 
-    sdata.add_shapes(SopaKeys.BAYSOR_BOUNDARIES, geo_df, overwrite=True)
+    sdata.shapes[SopaKeys.BAYSOR_BOUNDARIES] = geo_df
+    save_shapes(sdata, SopaKeys.BAYSOR_BOUNDARIES, overwrite=True)
 
     if sdata.table is not None:
         log.warn("Table already existing. It will be replaced by the new one.")

sopa/segmentation/patching.py

@@ -18,7 +18,13 @@
 from spatialdata.transformations import get_transformation
 
 from .._constants import EPS, ROI, SopaFiles, SopaKeys
-from .._sdata import get_boundaries, get_item, get_spatial_image, to_intrinsic
+from .._sdata import (
+    get_boundaries,
+    get_item,
+    get_spatial_image,
+    save_shapes,
+    to_intrinsic,
+)
 
 log = logging.getLogger(__name__)
 
@@ -204,7 +210,8 @@ def write(self, overwrite: bool = True, shapes_key: str | None = None) -> gpd.Ge
             geo_df, transformations=get_transformation(self.element, get_all=True)
         )
 
-        self.sdata.add_shapes(shapes_key, geo_df, overwrite=overwrite)
+        self.sdata.shapes[shapes_key] = geo_df
+        save_shapes(self.sdata, shapes_key, overwrite=overwrite)
 
         log.info(f"{len(geo_df)} patches were saved in sdata['{shapes_key}']")
 

sopa/segmentation/stainings.py

@@ -16,7 +16,7 @@
 from tqdm import tqdm
 
 from .._constants import SopaKeys
-from .._sdata import get_spatial_image
+from .._sdata import get_spatial_image, save_shapes
 from . import shapes
 
 log = logging.getLogger(__name__)
@@ -178,6 +178,7 @@ def add_shapes(cls, sdata: SpatialData, cells: list[Polygon], image_key: str, sh
         geo_df.index = image_key + geo_df.index.astype(str)
 
         geo_df = ShapesModel.parse(geo_df, transformations=get_transformation(image, get_all=True))
-        sdata.add_shapes(shapes_key, geo_df, overwrite=True)
+        sdata.shapes[shapes_key] = geo_df
+        save_shapes(sdata, shapes_key, overwrite=True)
 
         log.info(f"Added {len(geo_df)} cell boundaries in sdata['{shapes_key}']")

sopa/segmentation/tissue.py

@@ -10,7 +10,7 @@
 from spatialdata.models import ShapesModel
 
 from .._constants import ROI
-from .._sdata import get_intrinsic_cs, get_item
+from .._sdata import get_intrinsic_cs, get_item, save_shapes
 
 log = logging.getLogger(__name__)
 
@@ -113,6 +113,7 @@ def _save_tissue_segmentation(
         geo_df, image_scale.attrs["transform"][image_cs], maintain_positioning=True
     )
 
-    sdata.add_shapes(ROI.KEY, geo_df)
+    sdata.shapes[ROI.KEY] = geo_df
+    save_shapes(sdata, ROI.KEY)
 
     log.info(f"Tissue segmentation saved in sdata['{ROI.KEY}']")

sopa/utils/polygon_crop.py

@@ -13,7 +13,7 @@
 from spatialdata.transformations import get_transformation
 
 from .._constants import ROI
-from .._sdata import get_spatial_image
+from .._sdata import get_spatial_image, save_shapes
 from .image import resize
 
 log = logging.getLogger(__name__)
@@ -134,6 +134,7 @@ def polygon_selection(
     geo_df = ShapesModel.parse(
         geo_df, transformations=get_transformation(sdata[image_key], get_all=True)
     )
-    sdata.add_shapes(ROI.KEY, geo_df)
+    sdata.shapes[ROI.KEY] = geo_df
+    save_shapes(sdata, ROI.KEY)
 
     log.info(f"Polygon saved in sdata['{ROI.KEY}']")