---

title: "Conversion"

output: 

  html_document:

    toc: true

    toc_float:

      collapsed: false

      smooth_scroll: false

---

<style>

.title{

  display: none;

}

body {

  text-align: justify   

}

.center {

  display: block;

  margin-left: auto;

  margin-right: auto;

}

</style>

```{css, echo=FALSE}

.watch-out {

  color: black;

}

```

```{r setup, include=FALSE}

# use rmarkdown::render_site(envir = knitr::knit_global())

knitr::opts_chunk$set(highlight = TRUE, echo = TRUE)

```

<br>

# Conversion to Other Platforms

VoltRon is capable of end-to-end spatial data analysis for all levels of spatial resolutions, including those of single cell resolution. However, VoltRon provides a ecosystem friendly infrastructure where VoltRon objects could be transformed into data structures used by popular computational platforms such as [Seurat](https://satijalab.org/seurat/), [Squidpy](https://squidpy.readthedocs.io/en/stable/) and even [Zarr](http://vitessce.io/docs/data-file-types/#anndata-zarr) for interactive spatial data visualizatiob with [Vitessce](http://vitessce.io/).

For both **Seurat (R)** and **Squidpy (Python)**, we analyse readouts of the experiments conducted on example tissue sections analysed by the [Xenium In Situ](https://www.10xgenomics.com/platforms/xenium) platform. For more information on processing and analyzing Xenium datasets, check the [Cell/Spot Analysis](spotanalysis.html) tutorial.

<br>

## Seurat

We will first see how we can transform VoltRon objects into Seurat object and use built-in functions such as **FindAllMarkers** to visualize marker genes of clusters found by VoltRon. You can find the clustered Xenium data using VoltRon [here](https://bimsbstatic.mdc-berlin.de/landthaler/VoltRon/Cellanalysis/Xenium/Xenium_data_clustered.rds). 

```{r class.source="watch-out", eval = FALSE}

Xen_data <- readRDS("Xenium_data_clustered.rds")

SampleMetadata(Xen_data)

```

```

        Assay    Layer   Sample

Assay1 Xenium Section1 XeniumR1

Assay3 Xenium Section1 XeniumR2

```

<br>

We use the **as.Seurat** function to convert spatial assays of VoltRon into Seurat objects. Here, a Seurat object defines spatial components of cellular and subcellular assays as **FOV** objects, and we use the **type = "image"** argument to convert spatial coordinates of cells and molecules into individual FOV objects for each Xenium assay/layer in the VoltRon object.  

Please check the [Analysis of Image-based Spatial Data in Seurat](https://satijalab.org/seurat/articles/seurat5_spatial_vignette_2) tutorial for more information on analyzing FOV-based spatial data sets with Seurat. 

note: Use VoltRon::as.Seurat to avoid conflict with Seurat package's as.Seurat function

```{r class.source="watch-out", eval = FALSE}

library(Seurat)

Xen_data_seu <- VoltRon::as.Seurat(Xen_data, cell.assay = "Xenium", type = "image")

Xen_data_seu <- NormalizeData(Xen_data_seu)

Xen_data_seu

```

```

An object of class Seurat 

313 features across 283298 samples within 1 assay 

Active assay: Xenium (313 features, 0 variable features)

 1 layers present: counts

 2 dimensional reductions calculated: pca, umap

 2 spatial fields of view present: fov_Assay1 fov_Assay3

```

<br>

### Marker Analysis

Now that we converted VoltRon into a Seurat object, we can pick the **Clusters** metadata column indicating the clustering of cells and test for marker genes of each individual cluster. 

```{r class.source="watch-out", eval = FALSE}

Idents(Xen_data_seu) <- "Clusters"

markers <- FindAllMarkers(Xen_data_seu)

head(markers[order(markers$avg_log2FC, decreasing = TRUE),])

```

<div><pre><code style="font-size: 13px;">         p_val avg_log2FC pct.1 pct.2 p_val_adj cluster   gene

CPA3         0   7.343881 0.977 0.029         0      16   CPA3

CTSG         0   7.114698 0.878 0.011         0      16   CTSG

LILRA4.1     0   6.992717 0.939 0.015         0      19 LILRA4

ADIPOQ       0   6.860190 0.974 0.025         0       5 ADIPOQ

MS4A1        0   6.763083 0.919 0.027         0      17  MS4A1

BANK1        0   6.082192 0.889 0.037         0      17  BANK1</code></pre></div>

<br>

### Visualization

We can now pick top positive markers from each of these clusters prior to visualization.

```{r class.source="watch-out", eval = FALSE}

library(dplyr)

topmarkers <- markers %>%

  group_by(cluster) %>%

  top_n(n = 5, wt = avg_log2FC)

```

Here, VoltRon incorporates the unique markers learned by the **FindAllMarkers** function from Seurat and uses them to visualize the expression of these markers on heatmaps, and now we can also use these markers for annotating the clusters.   

```{r class.source="watch-out", eval = FALSE}

library(ComplexHeatmap)

marker_features <- unique(topmarkers$gene)

vrHeatmapPlot(Xen_data, features = marker_features, group.by = "Clusters", 

              show_row_names = TRUE, font.size = 10)

```

<img width="100%" height="100%" src="https://bimsbstatic.mdc-berlin.de/landthaler/VoltRon/Package/images/conversions_seurat_heatmap.png" class="center">

<br>

### Convert with Molecule Data

If defined, the **as.Seurat** function may also convert the molecule assay of the VoltRon object into a Seurat FOV object and allow visualizing molecules. You can find the Xenium VoltRon object with the molecule assay [here](https://bimsbstatic.mdc-berlin.de/landthaler/VoltRon/Cellanalysis/Xenium/Xen_R1.rds). 

```{r class.source="watch-out", eval = FALSE}

Xen_R1 <- readRDS("Xen_R1.rds")

SampleMetadata(Xen_R1)

```

```

            Assay    Layer   Sample

Assay1     Xenium Section1 XeniumR1

Assay2 Xenium_mol Section1 XeniumR1

```

<br>

We define both the cell level assay and the molecule level assay. 

```{r class.source="watch-out", eval = FALSE}

Xen_R1_seu <- as.Seurat(Xen_R1, cell.assay = "Xenium", molecule.assay = "Xenium_mol", type = "image")

```

<br>

Now we can visualize molecules alongside with cells.

```{r class.source="watch-out", eval = FALSE}

ImageDimPlot(Xen_R1_seu, fov = "fovAssay1", molecules = "PGR", group.by = "orig.ident", cols = "lightgrey", mols.size = 1)

```

<img width="60%" height="60%" src="https://bimsbstatic.mdc-berlin.de/landthaler/VoltRon/Package/images/conversions_seurat_imagedimplot.png" class="center">

<br>

## SpatialExperiment

VoltRon can also convert objects in [SpatialExperiment](https://www.bioconductor.org/packages/release/bioc/html/SpatialExperiment.html) objects. We are going to use the Xenium data clustered using VoltRon [here](https://bimsbstatic.mdc-berlin.de/landthaler/VoltRon/Cellanalysis/Xenium/Xenium_data_clustered.rds).

```{r class.source="watch-out", eval = FALSE}

Xen_data <- readRDS("Xenium_data_clustered.rds")

SampleMetadata(Xen_data)

```

We use the **as.SpatialExperiment** function to convert spatial assays of VoltRon into SpatialExperiment objects. Please check the [Introduction to the SpatialExperiment class](https://www.bioconductor.org/packages/release/bioc/vignettes/SpatialExperiment/inst/doc/SpatialExperiment.html) tutorial for more information. 

```{r class.source="watch-out", eval = FALSE}

library(SpatialExperiment)

spe <- as.SpatialExperiment(Xen_data, assay = "Xenium")

```

Here we can parse the image and visualize.

```{r class.source="watch-out", eval = FALSE}

img <- imgRaster(spe, 

                 sample_id = "Assay1", 

                 image_id = "main")

plot(img)

```

<br>

## Squidpy (Anndata, h5ad)

A true ecosystem friendly computational platform should support data types across multiple computing environments. By allowing users to convert VoltRon objects into annotated data matrix formats such as [anndata](https://github.com/scverse/anndata), we can use built-in spatial data analysis methods available on [squidpy](https://squidpy.readthedocs.io/en/stable/).

You can find the clustered and the annotated Xenium data using VoltRon [here](https://bimsbstatic.mdc-berlin.de/landthaler/VoltRon/Cellanalysis/Xenium/Xenium_data_clustered_annotated.rds). 

Anndata objects wrapped on h5ad files are commonly used by the [scverse](https://www.nature.com/articles/s41587-023-01733-8) ecosystem for single cell analysis which bring together numeruous tools maintained and distributed by a large community effort. Both squidpy and [scanpy](https://scanpy.readthedocs.io/en/stable/) are currently available on scverse. 

```{r class.source="watch-out", eval = FALSE}

Xen_data <- readRDS("Xenium_data_clustered_annotated.rds")

as.AnnData(Xen_data, assay = "Xenium", file = "Xen_adata_annotated.h5ad")

```

<br>

### Configure Squidpy (scverse)

Here, we use the [reticulate](https://rstudio.github.io/reticulate/) package to call **scverse** module in Python through a prebuilt anaconda environment. However, any python installation with the scverse module can be incorporated by reticulate.   

```{r class.source="watch-out", eval = FALSE}

library(reticulate)

use_condaenv("scverse", required = T)

```

We import some other necessary modules such as pandas, scanpy and squidpy.

```{python class.source="watch-out", eval = FALSE}

from pathlib import Path

import numpy as np

import pandas as pd

import matplotlib.pyplot as plt

import seaborn as sns

import scanpy as sc

import squidpy as sq

sc.logging.print_header()

```

<br>

### Filter and Normalize 

We read the annotated Xenium object that was saved as an h5ad file using the **as.Anndata** function in VoltRon, and process before analysis. For more information using scanpy and squidpy on Xenium datasets, check the [Analyzing Xenium data](https://squidpy.readthedocs.io/en/stable/notebooks/tutorials/tutorial_xenium.html) tutorial at squidpy webpage. 

```{python class.source="watch-out", eval = FALSE}

adata = sc.read_h5ad("Xen_adata_annotated.h5ad")

adata.layers["counts"] = adata.X.copy()

sc.pp.normalize_total(adata, inplace=True)

sc.pp.log1p(adata)

```

<br>

### Visualize

We use the **squidpy.pl.spatial_scatter** functions available in squidpy to visualize the spatial localization of cell types of both Xenium replicates.

```{python class.source="watch-out", eval = FALSE}

fig, ax = plt.subplots(1, 2, figsize=(10, 7))

sq.pl.spatial_scatter(adata, library_key = "library_id", library_id = "Assay1", 

                      color=["CellType"], shape=None, size=1, img = False, ax=ax[0])

sq.pl.spatial_scatter(adata, library_key = "library_id", library_id = "Assay3", 

                      color=["CellType"], shape=None, size=1, img = False, ax=ax[1])

plt.show(ax)

```

<img width="100%" height="100%" src="https://bimsbstatic.mdc-berlin.de/landthaler/VoltRon/Package/images/conversions_anndata_spatial_scatter.png" class="center">

<br>

### Neighborhood Enrichment

We can now use high level spatially-aware functions available in squidpy. We first establish spatial neighbors using the delaunay graphs. The spatial graph and distances will be stored under **.obsp** attribute/matrix. 

```{python class.source="watch-out", eval = FALSE}

sq.gr.spatial_neighbors(adata, coord_type="generic", delaunay=True)

adata

```

```

AnnData object with n_obs × n_vars = 283298 × 313

    obs: 'Count', 'Assay', 'Layer', 'Sample', 'Clusters', 'CellType', 'library_id'

    uns: 'log1p', 'spatial_neighbors'

    obsm: 'spatial'

    layers: 'counts'

    obsp: 'spatial_connectivities', 'spatial_distances'

```

We can now conduct the permutation test for neighborhood enrichment across cell type pairs. 

```{python class.source="watch-out", eval = FALSE}

sq.gr.nhood_enrichment(adata, cluster_key="CellType")

```

```{python class.source="watch-out", eval = FALSE}

fig, ax = plt.subplots(1, 2, figsize=(13, 7))

sq.pl.nhood_enrichment(adata, cluster_key="CellType", figsize=(8, 8), 

                       title="Neighborhood enrichment adata", ax=ax[0])

sq.pl.spatial_scatter(adata, color="CellType", library_key = "library_id", 

                      library_id = "Assay1", shape=None, size=2, ax=ax[1])

plt.show(ax)

```

<img width="100%" height="100%" src="https://bimsbstatic.mdc-berlin.de/landthaler/VoltRon/Package/images/conversions_anndata_neighborhood.png" class="center">

<br>

## Vitessce (Anndata, zarr) 

In this section, we will transform VoltRon objects of Xenium data into zarr arrays, and use them for interactive visualization in [Vitessce](http://vitessce.io/). We should first download the vitessceR package which incorporates wrapper function to visualize zarr arrays interactively in R.

```{r class.source="watch-out", eval = FALSE}

install.packages("devtools")

devtools::install_github("vitessce/vitessceR")

```

<br>

You can find the clustered and annotated Xenium data using VoltRon [here](https://bimsbstatic.mdc-berlin.de/landthaler/VoltRon/Cellanalysis/Xenium/Xenium_data_clustered_annotated.rds). 

```{r class.source="watch-out", eval = FALSE}

Xen_data <- readRDS("Xenium_data_clustered_annotated.rds")

SampleMetadata(Xen_data)

```

```

        Assay    Layer   Sample

Assay1 Xenium Section1 XeniumR1

Assay2 Xenium Section1 XeniumR2

```

<br>

### Interactive Visualization

Now we can convert the VoltRon object into a zarr array using the **as.Zarr** function which will create the array in a specified location.

```{r class.source="watch-out", eval = FALSE}

as.AnnData(Xen_data, assays = "Assay1", 

           file = "xendata_clustered_annotated.zarr", flip_coordinates = TRUE)

```

We can use the zarr file directly in the **vrSpatialPlot** function to visualize the zarr array interactively in Rstudio viewer. The **reduction** arguement allows the umap of the Xenium data to be visualized alongside with the spatial coordinates of the Xenium cells.

```{r class.source="watch-out", eval = FALSE}

vrSpatialPlot("xendata_clustered_annotated.zarr", group.by = "CellType", reduction = "umap")

```

<img width="100%" height="100%" src="https://bimsbstatic.mdc-berlin.de/landthaler/VoltRon/Package/images/conversions_interactive.png" class="center">

<br>

<img width="100%" height="100%" src="https://bimsbstatic.mdc-berlin.de/landthaler/VoltRon/Package/images/conversions_interactive_zoom.png" class="center">