How to extract transformed polygons from sdata.shapes?

Help spatialdata
1

I have some segmentations (of structures larger than cells) that were performed on an H&E image, that I have imported, added to my spatialdata object containing my Xenium data, and then applied the appropriate XeniumExplorer transformation matrix to.

I’d now like to find, for each polygon in the new segmentations, the IDs of the cells where the centroid falls inside the polygon.

I thought I’d be able to do this using the contains method of the Polygon object, but when I access the segmentation polygons using sdata.shapes['new_polygons'], the coordinates in the geometry column are the same as in my un-transformed input, and I don’t find any overlapping cells when I check the first couple of polygons.

Is it possible to extract/export the transformed coordinates of the polygons/points in sdata.shapes?

Or do you have another suggestion of how to do this?

I’ve looked at the polygon_query example, but I don’t follow what’s going on there - there’s clearly some kind of transformation happening between the polygon plotted by itself and the polygon plotted on top of the Visium data, but it’s not explained.

I tried that anyway, but I get an error, which I think is for the same reason as above - the untransformed polygon coords are being used for the query, so no cells overlap:

poly_cropped = polygon_query(
    sdata,
    polygon=sdata.shapes['tubule_polygons'].iloc[1,]["geometry"],
    target_coordinate_system="global",
)
File /mnt/storage/vaquerizas/liz/spatial/xenium/python_analysis/py3.10/lib/python3.10/site-packages/spatialdata/_core/query/relational_query.py:168, in _filter_table_by_elements(table, elements_dict, match_rows)
    166 assert set(elements_dict.keys()).issubset({"images", "labels", "shapes", "points"})
    167 assert len(elements_dict) > 0, "elements_dict must not be empty"
--> 168 assert any(
    169     len(elements) > 0 for elements in elements_dict.values()
    170 ), "elements_dict must contain at least one dict which contains at least one element"
    171 if table is None:
    172     return None

AssertionError: elements_dict must contain at least one dict which contains at least one element
2

Hi @liz-is, what you experiences is the expected behavior: we store the raw coordinates paired with lightweight transformation metadata and we keep the raw coordinates raw unless the user explicitly wants to transform them.

The solution in your case is to use transform() in both elements that you want to transform and then use the geopandas APIs (such as contains()). Please notice that you need to specify the same to_coordinate_system parameter for both calls, or leave the default value 'global'.

You may also find get_centroids() useful.

Two more hints:

  • if one of your segmentation data is raster, you can derive the polygons/multipolygons using to_polygons().

Finally, after you filter the cells you may want to also filter the table (which annotates the cells). For this we have the function match_table_to_element (for more control you can use join_spatialelement_table().

3

Thanks! I knew there had to be a simple way to do it.

I found the documentation for transform() slightly confusing, but I got it working! I think it’d be really helpful to have an example of using this function added to the tutorial on transformations at some point.

One remaining question - when I select and print a single polygon after transformation, its orientation doesn’t match how it’s plotted with sdata.pl.render_shapes() - is this expected?

4

Great to hear that you solved your problem, and thanks for the feedback on the documentation. I tracked it here.

Finally, the different orientation is expected as some plotting frameworks differ in how they orient the y axis; but you can quickly adjust this, for instance by plotting on a specific ax, and then using matplotlib APIs.