Updates tp day1-1b

Author: julien-roux

_freeze/day1/day1-1a_SpatialExperiment/execute-results/html.json

@@ -55,7 +55,8 @@ library(HDF5Array)
 
 We will work on a **Visium HD** dataset from a human colorectal cancer study by [Oliveira _et al_](https://www.nature.com/articles/s41588-025-02193-3). The dataset contains normal adjacent tissue (NAT) and colorectal carcinoma (CRC) from 5 patients.
 
-We will focus on a Visium HD slide from one patient (called "P2CRC"), available from the [10X website](https://www.10xgenomics.com/datasets/visium-hd-cytassist-gene-expression-libraries-of-human-crc). We will focus on the "binned output" available as an output of the Space Ranger pipeline (version 3.0). The full output is very large, so for practical reasons we use a region of interest, selected to include an interesting part of the tissue.:
+We will focus on a single Visium HD slide from a patient called "P2", available from the [10X website](https://www.10xgenomics.com/datasets/visium-hd-cytassist-gene-expression-libraries-of-human-crc). We will focus on the "binned output" available as an output of the Space Ranger pipeline (version 3.0). The full output is very large, so for practical reasons we use a region of interest, selected to include an interesting part of the tissue, because it contains gland-like epithelial structures together with stromal or lower-density areas:
+
 
 ```{r day1-1a-SpatialExperiment-4}
 roi_visium <- c(
@@ -106,7 +107,6 @@ if (!dir.exists("data/Human_Colon_Cancer_P2/")) {
 
 We will import the 16 um binned Visium HD output into a `SpatialExperiment` object and then subset it to the course region. The coordinates are in the native Visium HD image coordinate system returned by `spatialCoords()`.
 
-
 ```{r day1-1a-SpatialExperiment-6}
 spe <- TENxVisiumHD(
   spacerangerOut = "data/Human_Colon_Cancer_P2/",
@@ -210,7 +210,7 @@ rowData(spe) |> head()
 nrow(spe)
 ```
 
-The assay stores the count matrix:
+The assay stores the UMI counts matrix:
 
 ```{r day1-1a-SpatialExperiment-10}
 assay(spe)
@@ -246,7 +246,7 @@ imgRaster(spe) |> plot()
 
 ## Visualizing the tissue region
 
-To get a visual overview, we use `ggspavis::plotVisium`. First we plot the tissue image without bins:
+To get a visual overview, we use `plotVisium()` function from the `ggspavis` package. First we plot the tissue image without bins:
 
 ```{r day1-1a-SpatialExperiment-15}
 plotVisium(spe, spots = FALSE)
@@ -270,7 +270,7 @@ Why do we see only a specific region of the slide?
 The full Visium HD output is much larger than needed for this introductory exercise. At the start of this exercise, we subsetted the object to one representative tissue region so that it remains fast to plot and process. This same region is also used as the anchor for selecting a comparable Xenium region in `Exercise 1 B`.
 :::
 
-It is possible to colour bins according to a gene or a column in `colData`, for example `PIGR`:
+It is possible to colour bins according to a gene or a column in `colData`, for example `PIGR`. We renamed the rows to gene symbols after import, which is why we can refer to this marker as `"PIGR"` rather than by its Ensembl ID. At this stage we plot the number of UMIs since log-normalized counts will be introduced later in the course.
 
 ```{r day1-1a-SpatialExperiment-17}
 #| warning: false
@@ -339,7 +339,7 @@ We save the filtered `SpatialExperiment` object for the next steps. Remember tha
 ```{r day1-1a-SpatialExperiment-21}
 dir.create("results/day1", showWarnings = FALSE, recursive = TRUE)
 saveHDF5SummarizedExperiment(spe,
-  dir = "results/day1", prefix = "01.1_spe", replace = TRUE,
+  dir = "results/day1", prefix = "01.1_spe_", replace = TRUE,
   chunkdim = NULL, level = NULL, as.sparse = NA,
   verbose = NA
 )

day1/day1-1a_SpatialExperiment.qmd

@@ -30,17 +30,17 @@ library(Voyager)
 library(ggplot2)
 library(patchwork)
 library(scuttle)
+library(HDF5Array)
 ```
 
 ## Data for the course
 
-We will work with the **Xenium In Situ, sample P2 CRC** dataset from the same human colorectal cancer study used in `Exercise 1 A`. The Xenium section is a serial section of the Visium HD P2 CRC sample.
+We will work on a **Xenium** dataset from a human colorectal cancer study by [Oliveira _et al_](https://www.nature.com/articles/s41588-025-02193-3). The dataset contains normal adjacent tissue (NAT) and colorectal carcinoma (CRC) from 5 patients.
 
-This means the two datasets are biologically related but not the exact same physical tissue section. We therefore subset Xenium to an **approximately matching tissue region**, guided by the Visium HD region and by morphology. This makes the two exercises comparable while avoiding a full image-registration workflow during the practical.
+We will focus on a single Xenium slide from a patient called "P2", available from the [10X website](https://www.10xgenomics.com/platforms/visium/product-family/dataset-human-crc). This slide is a serial section of the sample used in `Exercise 1 A`, so the two datasets are biologically related but not the exact same physical tissue section. 
 
-The full Xenium output is large and includes cell boundaries, nucleus boundaries, transcript locations, morphology images, and expression matrices. For this exercise, we use a compact course data folder that contains the required Xenium output files.
+The full Xenium output is very large and includes cell boundaries, nucleus boundaries, transcript locations, morphology images, and expression matrices. For this exercise, for practical reasons we use a region of interest, selected to **approximately match the region of interest used in `Exercise 1 A`**, and only include some of the layers in the imported compact course folder. This makes the two exercises comparable while avoiding a full image-registration workflow during the practical. This region was chosen from the overview plot because it contains gland-like epithelial structures together with stromal or lower-density areas.
 
-The selected Xenium region is:
 
 ```{r day1-1b-SFE-2}
 roi_xenium <- c(
@@ -52,9 +52,7 @@ roi_xenium <- c(
 roi_xenium
 ```
 
-We import the compact Xenium output with `flip = "none"` so the exercise uses the native Xenium coordinate orientation. The same visual region may appear with negative y coordinates if the geometry has been flipped during import with morphology images. Because the compact course folder does not include morphology images, we use the positive y coordinates shown above.
-
-This region was chosen from the overview plot because it contains gland-like epithelial structures together with stromal or lower-density areas.
+We import the Xenium output with `flip = "none"` so the exercise uses the native Xenium coordinate orientation. The same visual region may appear with negative y coordinates if the geometry has been flipped during import with morphology images. Because the compact course folder does not include morphology images, we use the positive y coordinates shown above.
 
 ```{r day1-1b-SFE-3}
 options(timeout = 600)
@@ -83,10 +81,6 @@ We import the compact Xenium output into a `SpatialFeatureExperiment` object:
 ```{r day1-1b-SFE-4}
 sfe_full <- readXenium("data/Human_Colon_Cancer_P2/xenium/outs", flip = "none")
 
-if (ncol(sfe_full) == 0) {
-  stop("The Xenium object has 0 cells before ROI subsetting. Check that cell_feature_matrix.h5 is present in data/Human_Colon_Cancer_P2/xenium/outs/.")
-}
-
 rownames(sfe_full) <- uniquifyFeatureNames(
   ID = rowData(sfe_full)$ID,
   names = rowData(sfe_full)$Symbol
@@ -97,27 +91,23 @@ range(spatialCoords(sfe_full)[, 1])
 range(spatialCoords(sfe_full)[, 2])
 ```
 
-Then we subset it to the course region:
+Then we subset it to the region of interest:
 
 ```{r day1-1b-SFE-5}
 sfe <- sfe_full[, spatialCoords(sfe_full)[, 1] >= roi_xenium[["xmin"]] &
   spatialCoords(sfe_full)[, 1] <= roi_xenium[["xmax"]] &
   spatialCoords(sfe_full)[, 2] >= roi_xenium[["ymin"]] &
   spatialCoords(sfe_full)[, 2] <= roi_xenium[["ymax"]]]
 
-if (ncol(sfe) == 0) {
-  stop("The Xenium ROI selected 0 cells. Check range(spatialCoords(sfe)) and revise roi_xenium.")
-}
-
 sfe
 ```
 
 ::: callout-important
 ## Exercise 1
 
-Why is a `SpatialFeatureExperiment` useful for Xenium data?
+Why is an object of the `SpatialFeatureExperiment` class used to store Xenium data?
 
-Which information does Xenium provide that is not naturally represented by a simple spot- or bin-level `SpatialExperiment`?
+Which information could not be represented by a simple spot- or bin-level `SpatialExperiment`?
 :::
 
 ::: {.callout-tip collapse="true"}
@@ -130,7 +120,7 @@ Xenium is an image-based spatial transcriptomics technology. In addition to a co
 
 ## Exploring the object
 
-`SpatialFeatureExperiment` extends `SpatialExperiment`, so familiar accessors still apply:
+`SpatialFeatureExperiment` extends `SpatialExperiment`, which extends `SingleCellExperiment`, so familiar accessors still apply:
 
 ```{r day1-1b-SFE-6}
 dim(sfe)
@@ -148,24 +138,33 @@ Questions:
 
 -   What do the columns represent?
 -   What do the rows represent?
--   How is this different from the Visium HD `spe` object in Exercise 1 A?
+-   What differs compared to the Visium HD `spe` object in `Exercise 1 A`?
+-   What is stored in the main assay?
 :::
 
 ::: {.callout-tip collapse="true"}
 ## Answer
 
-In this object, columns represent Xenium cells, while rows represent measured genes/features. This differs from the Visium HD binned `SpatialExperiment`, where columns represent spatial bins rather than segmented cells.
+In this object, columns represent segmented cells, while rows represent genes. This differs from the Visium HD slide stored in a `SpatialExperiment` object, where columns represent spatial bins rather than segmented cells.
 
 ```{r day1-1b-SFE-7}
-dim(sfe)
 ncol(sfe)
 nrow(sfe)
 ```
+
+We notice that the number of features is small (541) as the Xenium technology measures a targeted panel rather than the whole transcriptome.
+
+```{r day1-1b-SFE-7b}
+assay(sfe)
+```
+
+The count matrix indicates the number of observed molecules for each of the 541 genes in each cell. This matrix is analogous to a UMI counts matrix in the Visium HD dataset.
+
 :::
 
 ## Inspecting geometries
 
-One of the key features of `SpatialFeatureExperiment` is that it can store geometries associated with cells or features.
+In addition to the slots from the `SpatialExperiment` object, one of the key features of `SpatialFeatureExperiment` is that it can store geometries associated with cells or features.
 
 ```{r day1-1b-SFE-8}
 colGeometries(sfe)
@@ -174,24 +173,24 @@ colGeometries(sfe)
 ::: callout-important
 ## Exercise 3
 
-Inspect the available column geometries in the Xenium object.
-
-Which geometries are available? What biological structures do they correspond to?
+Inspect the available column geometries in the Xenium object. To which biological structures do they correspond to? In which class are stored these objects?
 :::
 
 ::: {.callout-tip collapse="true"}
 ## Answer
 
-For Xenium data, the available geometries commonly include cell centroids, cell boundaries, and nucleus boundaries. The exact names depend on the reader and package versions.
+The available geometries include cell centroids, cell boundaries, and nucleus boundaries. These are are commonly available for imaging-based spatial transcriptomics technologies, but the exact names depend on the reader functions and packages used to import the data.
 
 ```{r day1-1b-SFE-9}
-colGeometries(sfe)
+colGeometries(sfe)[["centroids"]]
 ```
+The geometries are stored using "Simple Features" objects, which allow to encode spatial data and can be easily manipulated with the [`sf` package](https://r-spatial.github.io/sf/index.html)
+
 :::
 
 ## Visualizing Xenium cells
 
-We can start with a simple coordinate plot of the selected Xenium region:
+We can start with a simple scatterplot of the cells in the selected region of interest:
 
 ```{r day1-1b-SFE-10}
 xy <- as.data.frame(spatialCoords(sfe))
@@ -207,7 +206,7 @@ ggplot(xy, aes(x = x, y = y)) +
 ::: callout-important
 ## Exercise 4
 
-Compare this point-based visualization with the Visium HD plot from Exercise 1 A.
+Compare this point-based visualization with the Visium HD plot from `Exercise 1 A`.
 
 What does one point represent in each dataset?
 :::
@@ -218,9 +217,7 @@ What does one point represent in each dataset?
 In the Visium HD exercise, each point or square represents a 16 um spatial bin. In the Xenium exercise, each point represents a segmented cell centroid. This is one of the key conceptual differences between binned sequencing-based spatial data and image-based cell-resolved spatial data.
 :::
 
-We can also visualize expression for a marker gene. For example, `PIGR` is present in the Xenium panel and shows a clear spatial pattern in this region. We renamed the rows to gene symbols after import, which is why we can refer to this marker as `"PIGR"` rather than by its Ensembl ID.
-
-At this stage of the course we plot `raw counts`. `Log-normalized` values will be introduced later, in the normalization part of the course.
+We can also visualize expression for a marker gene. For example, `PIGR` is present in the Xenium panel and shows a clear spatial pattern in this region. We renamed the rows to gene symbols after import, which is why we can refer to this marker as `"PIGR"` rather than by its Ensembl ID. At this stage we plot the number of molecules since log-normalized counts will be introduced later in the course. The `plotSpatialFeature()` function from the `Voyager` package is handy:
 
 ```{r day1-1b-SFE-11}
 plotSpatialFeature(sfe, "PIGR", exprs_values = "counts")
@@ -240,47 +237,34 @@ How does the interpretation differ from plotting a gene on Visium HD bins?
 The Xenium signal is measured over segmented cells, while the Visium HD binned signal aggregates molecules within fixed spatial bins. Xenium can therefore represent cell-level heterogeneity more directly, but it measures a targeted panel rather than the whole transcriptome.
 :::
 
-## Subsetting with coordinates
-
-At the start of the exercise, we subsetted the full imported Xenium object to a course region. We can still apply coordinate-based subsetting again in the same spirit:
-
-```{r day1-1b-SFE-12}
-xy <- spatialCoords(sfe)
-
-roi <- c(
-  xmin = unname(quantile(xy[, 1], 0.25)),
-  xmax = unname(quantile(xy[, 1], 0.75)),
-  ymin = unname(quantile(xy[, 2], 0.25)),
-  ymax = unname(quantile(xy[, 2], 0.75))
-)
-
-keep <- xy[, 1] > roi[["xmin"]] & xy[, 1] < roi[["xmax"]] &
-  xy[, 2] > roi[["ymin"]] & xy[, 2] < roi[["ymax"]]
-
-sfe_sub <- sfe[, keep]
-dim(sfe)
-dim(sfe_sub)
-```
-
-::: callout-important
 ## Exercise 6
 
-Why did we subset Xenium at the start of the exercise instead of using the full public output directly throughout the practical?
-:::
+Use the `plotSpatialFeature()` function to visualize the cell segmentation mask, colored by cell area, and the nuclei segmentation mask, colored by nucleus area
 
 ::: {.callout-tip collapse="true"}
 ## Answer
 
-The full Xenium output is large and contains many cells, boundaries, transcript coordinates, and morphology images. A course-sized subset keeps the exercise responsive while preserving the structure of a real Xenium object.
+```{r day1-1b-SFE-12}
+plotSpatialFeature(sfe, 
+                   colGeometryName="cellSeg", 
+                   features="cell_area") 
 
-The subset is also chosen to approximately match the Visium HD P2 CRC region from Exercise 1 A. Because Visium HD and Xenium are serial sections, this is a morphology-matched comparison rather than an exact same-cell comparison.
+plotSpatialFeature(sfe, 
+                   colGeometryName="nucSeg", 
+                   features="nucleus_area") 
+```
 :::
 
 ## Save the object
 
 ```{r day1-1b-SFE-13}
 dir.create("results/day1", showWarnings = FALSE, recursive = TRUE)
-saveRDS(sfe, "results/day1/01.1b_sfe_p2.rds")
+
+saveHDF5SummarizedExperiment(sfe,
+  dir = "results/day1", prefix = "01.1b_sfe_", replace = TRUE,
+  chunkdim = NULL, level = NULL, as.sparse = NA,
+  verbose = NA
+)
 ```
 
 Clear your environment: