Skip to content

Q: How to render xenium morphology image with original colors? #370

Description

@tsvvas

Hello @LucaMarconato!

Thank you very much for this tool. I want to plot a part of morphology image acquired with xenium + multimodal segmentation kit. It is loaded into spatialdata and plotted without errors, but the colors are dim and they are different from what I'd expected (DAPI - #0F73E6, Boundary - #F300A5, Interior RNA - #A4A400, interior protein - #008A00, reference values from Xenium Explorer). I tried passing a list of colors like this:

crop = lambda sdata: spatialdata.bounding_box_query(
    sdata,
    min_coordinate=[17_500, 55_000],
    max_coordinate=[19_500, 57_000],
    axes=("x", "y"),
    target_coordinate_system="global",
)

crop(sdata).pl.render_images("morphology_focus", ["#0F73E6", "#F300A5", "#A4A400", "#008A00"]).pl.show(
    ax=axes[0], title="Morphology image", coordinate_systems="global"
)

The code does not produce an error, it works twice as long, and it results in the same image with colors unchanged.

The key question is how to adjust the colors and saturation of the plot, so that it looked closer to what I see in Xenium Explorer?

Secondary question is how to translate physical coordinates into spatialdata's global coordinate system? Am I supposed to get transformation and then apply transformation to a single or a couple points to get min and max for cropping above? I can make a separate issue with this question, if necessary.

Best,
Vasily

Activity

  1. LucaMarconato commented on Oct 10, 2024

    @LucaMarconato
    Member

    Hi, for coloring the channels please use the palette argument.

    Regarding coordinates yes, once you identify the transformation from the physical coordinate system to the global one, you would need to transform the coordinates and use these in the bounding box query.

    We are going to improve the ergonomics of this as described in these issues (except for 338, which is unrelated): https://github.com/scverse/spatialdata/issues?q=sort%3Aupdated-desc%20is%3Aissue%20is%3Aopen%20physical

  2. LucaMarconato commented on Oct 10, 2024

    @LucaMarconato
    Member

    Also please check my comment here (which answers almost a duplicate of question 2): scverse/spatialdata-io#205

  3. timtreis commented on Oct 15, 2024

    @timtreis
    Member

    Hey @tsvvas,

    I assume you have an image with 4 channels that you're trying to map the individual colors to, right? As Luca said, you can use the palette argument for this, see here:

    Image

    Please let us know if that solves your issue :)

  4. tsvvas commented on Oct 15, 2024

    @tsvvas
    Author

    Hi @timtreis,

    Thank you for checkin in on the issue. I managed to produce an image, but not without issues.

    First, the morphology_focus image after reading with spatiadata_io.xenium() has 5 channels for some reason, and I don't know what is the 5th channel. To my best knowledge there should be only 4 images.

    Second, there is a small issue where providing only 4 colors causes an error, even if I select only 4 channels. So to produce the image I had to pass five colors like this:

    cropped.pl.render_images(
        "morphology_focus",
        palette=["#0F73E6", "#F300A5", "#A4A400", "#008A00", "#FFFFFF"],
        channels=[0, 1, 2, 3],
    ).pl.show(
        title="Morphology image", 
        coordinate_systems="global",
    )

    At first I thought that the rendered image is too dim, but it appears that one of the layers is not rendered properly this way. Here is the image I get after executing the code above:

    Here is the image of the same region from Xenium Explorer:

    And here is the image of the same region from Xenium Explorer without Interior RNA:

    So, I am still looking for the way to render the image the same way it is shown in Xenium Explorer.

    Best,
    Vasily

  5. timtreis commented on Oct 15, 2024

    @timtreis
    Member

    The argument is actually called channel without the s, could you try that? It requires 5 colors since you're not limiting the channels and therefore it uses all 5 you have in the image.

  6. tsvvas commented on Oct 15, 2024

    @tsvvas
    Author

    With channel I get the following error:

    cropped.pl.render_images(
        "morphology_focus",
        palette=["#0F73E6", "#F300A5", "#A4A400", "#008A00"],
        channel=[0, 1, 2, 3],
    ).pl.show(
        title="Morphology image", 
        coordinate_systems="global"
    )

    Error:

    File [spatialdata_plot/pl/basic.py:862]
        856     if wanted_images_on_this_cs:
        857         rasterize = (params_copy.scale is None) or (
        858             isinstance(params_copy.scale, str)
        859             and params_copy.scale != "full"
        860             and (dpi is not None or figsize is not None)
        861         )
    --> 862         _render_images(
        863             sdata=sdata,
        864             render_params=params_copy,
        865             coordinate_system=cs,
        866             ax=ax,
        867             fig_params=fig_params,
        868             scalebar_params=scalebar_params,
        869             legend_params=legend_params,
        870             rasterize=rasterize,
        871         )
    File [spatialdata_plot/pl/render.py:634]
        632 layers = {}
        633 for ch_index, c in enumerate(channels):
    --> 634     layers[c] = img.sel(c=c).copy(deep=True).squeeze()
        636     if not isinstance(render_params.cmap_params, list):
        637         if render_params.cmap_params.norm is not None:
    File [xarray/core/dataarray.py:1670]
       1554 def sel(
       1555     self,
       1556     indexers: Mapping[Any, Any] | None = None,
       (...)
       1560     **indexers_kwargs: Any,
       1561 ) -> Self:
       1562     """Return a new DataArray whose data is given by selecting index
       1563     labels along the specified dimension(s).
       1564 
       (...)
       1668     Dimensions without coordinates: points
       1669     """
    -> 1670     ds = self._to_temp_dataset().sel(
       1671         indexers=indexers,
       1672         drop=drop,
       1673         method=method,
       1674         tolerance=tolerance,
       1675         **indexers_kwargs,
       1676     )
       1677     return self._from_temp_dataset(ds)
    File [xarray/core/dataset.py:3184]
       3116 """Returns a new dataset with each array indexed by tick labels
       3117 along the specified dimension(s).
       3118 
       (...)
       3181 
       3182 """
       3183 indexers = either_dict_or_kwargs(indexers, indexers_kwargs, "sel")
    -> 3184 query_results = map_index_queries(
       3185     self, indexers=indexers, method=method, tolerance=tolerance
       3186 )
       3188 if drop:
       3189     no_scalar_variables = {}
    File [xarray/core/indexing.py:193]
        191         results.append(IndexSelResult(labels))
        192     else:
    --> 193         results.append(index.sel(labels, **options))
        195 merged = merge_sel_results(results)
        197 # drop dimension coordinates found in dimension indexers
        198 # (also drop multi-index if any)
        199 # (.sel() already ensures alignment)
    File [xarray/core/indexes.py:791]
        789                 indexer = self.index.get_loc(label_value)
        790             except KeyError as e:
    --> 791                 raise KeyError(
        792                     f"not all values found in index {coord_name!r}. "
        793                     "Try setting the `method` keyword argument (example: method='nearest')."
        794                 ) from e
        796 elif label_array.dtype.kind == "b":
        797     indexer = label_array
    
    KeyError: "not all values found in index 'c'. Try setting the `method` keyword argument (example: method='nearest')."
  7. tsvvas commented on Oct 15, 2024

    @tsvvas
    Author

    A fully reproducible example:

    ! curl -O https://cf.10xgenomics.com/samples/xenium/2.0.0/Xenium_V1_humanLung_Cancer_FFPE/Xenium_V1_humanLung_Cancer_FFPE_outs.zip --output-dir ../data
    ! mkdir ../data/Xenium_V1_humanLung_Cancer_FFPE_outs
    ! unzip ../data/Xenium_V1_humanLung_Cancer_FFPE_outs.zip -d ../data/Xenium_V1_humanLung_Cancer_FFPE_outs
    
    DATA_PATH = "../data/Xenium_V1_humanLung_Cancer_FFPE_outs"
    sdata = spio.xenium(DATA_PATH, n_jobs=-1)
    
    OUT_PATH = f"../output/{date.today()}"
    FILE_PATH = os.path.join(OUT_PATH, "data.zarr")
    os.mkdir(OUT_PATH)
    sdata.write(FILE_PATH)
    
    sdata = spatialdata.read_zarr(FILE_PATH)
    
    sdata.pl.render_images(
        "morphology_focus",
        palette=["#0F73E6", "#F300A5", "#A4A400", "#008A00"],
        channel=[0, 1, 2, 3],
    ).pl.show(
        title="Morphology image", 
        coordinate_systems="global"
    )

    Interestingly, if I provide five colors, the error becomes:
    AttributeError: 'DataTree' object has no attribute 'c'

    >>> spatialdata.__version__
    '0.2.3'
  8. timtreis commented on Oct 15, 2024

    @timtreis
    Member

    Hm, I think the true underlying issue here is that the image is stored as a datatree (multi-scale image) and we're not respecting that when sub-selecting for channels. Will check your MRE later, thank you!

    If you want, you could try to overwrite the multi-scale image with one of the scales of the datatree. It probably looks like the second one here, right? Image

  9. tsvvas commented on Oct 15, 2024

    @tsvvas
    Author

    Yes, this is how it looks like in the example above

    Image

    What is the 5th channel here? Or should I ask this in the spatialdata_io issues?

  10. timtreis commented on Oct 16, 2024

    @timtreis
    Member

    I'll probably create a new issue for the issue you probably identified, that we cannot pick channels out of datatrees robustly. So we might close this afterwards.

    Could you plot the individual channels of your image separately before opening an issue in spatialdata-io? I'd be surprised if we create a non-sense channel during reading. The separate plot of the 5 channels might help troubleshoot here

  11. tsvvas commented on Oct 18, 2024

    @tsvvas
    Author

    Hi @timtreis,

    Seems like the 5th channel is empty:

    >>> sdata.images["morphology_focus"]["scale0"]["image"].isel(c=4).sum(dim=["y", "x"]).values
    array(0, dtype=uint64)
  12. tsvvas commented on Oct 18, 2024

    @tsvvas
    Author

    I tried the overwriting with single-scale image approach and still got the same error (continuing the example with the data above):

    img = sdata.images["morphology_focus"].copy()
    img = img.drop_nodes(["scale1", "scale2", "scale3", "scale4"])
    sdata.images["morpology_single_scale"] = img
    sdata.pl.render_images(
        "morpology_single_scale",
        palette=["#0F73E6", "#F300A5", "#A4A400", "#008A00"],
        channels=[0, 1, 2, 3],
    ).pl.show(
        title="Morphology image", 
        coordinate_systems="global"
    )

    The error is still ValueError: If 'palette' is provided, its length must match the number of channels.

    Providing 5 colors in the palette still results in the absence of Interior RNA channel in the rendered image.

  13. timtreis commented on Oct 19, 2024

    @timtreis
    Member

    @LucaMarconato have we seen a 5th channel popping up before?

  14. tsvvas commented on Oct 22, 2024

    @tsvvas
    Author

    @timtreis The data from spatialdata tutorial on xenium analysis has 5 channels in "morphology_focus" image. But for me the key problem now is that I cannot render properly the other 4 channels.

  15. LucaMarconato commented on Oct 23, 2024

    @LucaMarconato
    Member

    @LucaMarconato have we seen a 5th channel popping up before?

    Yes, adding a dummy fifth channel was a workaround when the support for the rgb parameter in spatial-image (and napari-spatialdata) was not available. We can remove that from the xenium() reader now, I'll track this in a separate issue.

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment

Metadata

Metadata

Assignees

No one assigned

    Type

    No type

    Projects

    No projects

      Milestone

      No milestone

      Relationships

      None yet

      Development

      No branches or pull requests

      Issue actions