Version cleanup

conda-scimap.yaml

@@ -0,0 +1,43 @@
+# Setup:
+#   mamba env create -f conda-scimap.yaml
+#   mamba activate spatialcells-scimap
+#   pip install --no-deps .
+#   pip install --no-deps scimap==2.3.6 smart_open
+
+# scimap pins numba 0.58 / llvmlite 0.41 / scipy 1.12 / pandas <2.2 / numpy <2,
+# which strictly forces the rest of the stack to certain older versions captured here.
+
+name: spatialcells-scimap
+
+channels:
+  - conda-forge
+
+dependencies:
+  - python >=3.10,<3.12
+  - pip
+  - ipykernel
+  - nbconvert
+  - nbclient
+
+  # Core runtime (pinned to scimap-compatible versions)
+  - shapely >=2.0
+  - numpy <2
+  - pandas <2.2
+  - scipy =1.12.*
+  - scikit-learn >=1.6
+  - matplotlib <3.9
+  - seaborn
+  - tqdm
+  - leidenalg
+  - scanpy <1.11
+  - anndata <0.11
+
+  # Optional: for joint analyses with other packages (as was done in the tutorials)
+  - squidpy =1.4.*
+  - plotly
+  - numba =0.58.*
+  - llvmlite =0.41.*
+  - dask >=2023.11,<2024
+  - dask-image <2024
+  - zarr >=2.11,<3
+  - setuptools <81

conda.yaml

@@ -1,19 +1,35 @@
-name: spatial-cells-env
+# Setup:
+#   mamba env create -f conda.yaml
+#   mamba activate spatialcells
+#   pip install --no-deps .
+
+# Use conda-scimap.yaml instead if you need scimap in the same env as well.
+
+name: spatialcells
 
 channels:
   - conda-forge
-  - defaults
 
 dependencies:
-  - python>=3.7, <3.11
+  - python >=3.10,<3.13
+  - pip
   - ipykernel
-  - matplotlib>=3.7.0
-  - pandas>=2.0.3
-  - seaborn>=0.12.2
-  - shapely>=2.0
+  - nbconvert
+  - nbclient
+
+  # Core runtime
+  - shapely >=2.0 
+  - numpy
+  - pandas
+  - scipy
+  - scikit-learn
+  - matplotlib
+  - seaborn
   - tqdm
-  - pip
-  - pip:
-      - anndata==0.9.2
-      - scanpy==1.9.4
-
+  - leidenalg
+  - scanpy
+  - anndata
+
+  # Optional: for joint analyses with other packages (as was done in the tutorials)
+  - squidpy
+  - plotly

spatialcells/__init__.py

@@ -1,8 +1,8 @@
-from pkg_resources import get_distribution, DistributionNotFound
+from importlib.metadata import version, PackageNotFoundError
 
 try:
-    __version__ = get_distribution("spatialcells").version
-except DistributionNotFound:
+    __version__ = version("spatialcells")
+except PackageNotFoundError:
     __version__ = "(local)"
 
 

spatialcells/measurements/_getSlidingWindowsComposition.py

@@ -1,3 +1,4 @@
+import anndata as ad
 import numpy as np
 import pandas as pd
 
@@ -26,33 +27,54 @@ def getSlidingWindowsComposition(
     :returns: A dataframe containing the cell type composition of the region in each window
     """
     if region_subset is None:
-        cells_roi = adata
+        obs = adata.obs
     else:
-        cells_roi = adata[adata.obs[region_col].isin(region_subset)]
-    cells_roi_maxx = int(cells_roi.obs["X_centroid"].max())
-    cells_roi_maxy = int(cells_roi.obs["Y_centroid"].max())
-    cells_roi_minx = int(cells_roi.obs["X_centroid"].min())
-    cells_roi_miny = int(cells_roi.obs["Y_centroid"].min())
-    all_windows_comp_df = []
-    for x in range(cells_roi_minx, cells_roi_maxx + window_size, step_size):
-        for y in range(cells_roi_miny, cells_roi_maxy + window_size, step_size):
-            cells_roi_window = cells_roi[
-                (cells_roi.obs["X_centroid"] >= x)
-                & (cells_roi.obs["X_centroid"] < x + window_size)
-                & (cells_roi.obs["Y_centroid"] >= y)
-                & (cells_roi.obs["Y_centroid"] < y + window_size)
-            ]
-            if cells_roi_window.shape[0] > min_cells:
-                cells_roi_window_composition = getRegionComposition(
-                    cells_roi_window, phenotype_col
-                )
-                cells_roi_window_composition["X_start"] = x
-                cells_roi_window_composition["Y_start"] = y
-                cells_roi_window_composition["window_size"] = window_size
-                cells_roi_window_composition["step_size"] = step_size
-                all_windows_comp_df.append(cells_roi_window_composition)
-    all_windows_comp_df = pd.concat(all_windows_comp_df)
-    return all_windows_comp_df
+        obs = adata.obs[adata.obs[region_col].isin(region_subset)]
+    if obs.empty:
+        return pd.DataFrame()
+
+    minx = int(obs["X_centroid"].min())
+    miny = int(obs["Y_centroid"].min())
+    maxx = int(obs["X_centroid"].max())
+    maxy = int(obs["Y_centroid"].max())
+
+    def composition_for(window_obs, x, y):
+        comp = getRegionComposition(
+            ad.AnnData(obs=window_obs), phenotype_col, regioncol=region_col
+        )
+        comp["X_start"] = x
+        comp["Y_start"] = y
+        return comp
+
+    rows = []
+    if window_size == step_size:
+        # speed up: each cell falls into exactly one window, so bucket once
+        # and group, rather than iterating all (x, y) pairs
+        x_bin = ((obs["X_centroid"] - minx) // step_size * step_size + minx).astype(int)
+        y_bin = ((obs["Y_centroid"] - miny) // step_size * step_size + miny).astype(int)
+        for (x, y), window_obs in obs.assign(_X_bin=x_bin, _Y_bin=y_bin).groupby(
+            ["_X_bin", "_Y_bin"], observed=True
+        ):
+            if window_obs.shape[0] > min_cells:
+                rows.append(composition_for(window_obs, int(x), int(y)))
+    else:
+        for x in range(minx, maxx + window_size, step_size):
+            for y in range(miny, maxy + window_size, step_size):
+                window_obs = obs[
+                    (obs["X_centroid"] >= x)
+                    & (obs["X_centroid"] < x + window_size)
+                    & (obs["Y_centroid"] >= y)
+                    & (obs["Y_centroid"] < y + window_size)
+                ]
+                if window_obs.shape[0] > min_cells:
+                    rows.append(composition_for(window_obs, x, y))
+
+    if not rows:
+        return pd.DataFrame()
+    out = pd.concat(rows)
+    out["window_size"] = window_size
+    out["step_size"] = step_size
+    return out
 
 
 def get_comp_mask(df, pheno_col, pheno_vals, step_size):

spatialcells/spatial/_getBoundary.py

@@ -1,3 +1,4 @@
+import numpy as np
 import shapely
 from shapely.geometry import Polygon, MultiPolygon
 from shapely.validation import make_valid
@@ -29,26 +30,31 @@ def getBoundary(anndata, communitycolumn, communityIndexList, alpha=100, debug=F
     xy = anndata.obs[anndata.obs[communitycolumn].isin(communityIndexList)][
         ["X_centroid", "Y_centroid"]
     ].to_numpy()
+    # slightly nudge the input points so they lie strictly inside.
+    eps = max(float(np.ptp(xy, axis=0).max()), 1.0) * 1e-9
     edge_components = getAlphaShapes(xy, alpha)
     polygons = [Polygon(edge[:, 0, :]).buffer(0) for edge in edge_components]
     boundary = MultiPolygon(polygons)
-    # make sure the boundary is valid and points of interest are inside
-    boundary = make_valid(boundary).buffer(1)
+    boundary = make_valid(boundary).buffer(eps)
     if boundary.geom_type == "Polygon":
         boundary = MultiPolygon([boundary])
+
     new_boundary = []
     for poly in boundary.geoms:
         holes = poly.interiors
         hole_polygons = [Polygon(hole) for hole in holes]
-        # prevent holes from touching the boundary
-        hole_limit = Polygon(poly.exterior).buffer(-1)
+        # shrink the limit by eps avoid the edge case of holes (internal polygons)
+        # being exactly tangential to the exterior, leading to an invalid polygon. 
+        hole_limit = Polygon(poly.exterior).buffer(-eps)
         hole_multipoly = shapely.unary_union(hole_polygons) & hole_limit
         if hole_multipoly.geom_type == "Polygon":
             hole_multipoly = MultiPolygon([hole_multipoly])
         new_holes_polygons = [p.exterior.coords for p in hole_multipoly.geoms]
         new_poly = Polygon(poly.exterior, new_holes_polygons)
         new_boundary.append(new_poly)
     boundary = make_valid(MultiPolygon(new_boundary))
+    if boundary.geom_type == "Polygon":
+        boundary = MultiPolygon([boundary])
 
     if debug:
         return boundary, polygons, edge_components