diff --git a/02_session2.Rmd b/02_session2.Rmd index afcdba7..8e9428d 100644 --- a/02_session2.Rmd +++ b/02_session2.Rmd @@ -723,6 +723,9 @@ spatInSituPlotPoints(visiumHD, knitr::include_graphics("img/02_session2/21-spatInSituPlotPoints.png") ``` + + + ## Hexbin 25 Goal is to create a higher resolution bin (hex25) and add to the Giotto object. We will aim to identify individual cell types and local neighborhood niches. @@ -740,6 +743,50 @@ visiumHD_subset = subsetGiottoLocs(gobject = visiumHD, y_max = 45500) ``` + +```{r, eval=FALSE} +spatInSituPlotPoints(visiumHD_subset, + show_image = F, + feats = NULL, + show_legend = F, + spat_unit = 'hex100', + point_size = 0.5, + show_polygon = TRUE, + use_overlap = FALSE, + polygon_feat_type = 'hex100', + polygon_fill_as_factor = TRUE, + polygon_fill = 'leiden_clus', + polygon_color = 'black', + polygon_line_size = 0.3) +``` + +```{r, echo=FALSE, out.width="100%", fig.align="center", fig.cap="xxx"} +knitr::include_graphics("img/02_session2/22-spatInSituPlotPoints.png") +``` + +```{r, eval=FALSE} +spatInSituPlotPoints(visiumHD_subset, + show_image = F, + feats = list('rna' = selected_feats), + feats_color_code = my_colors, + show_legend = F, + spat_unit = 'hex100', + point_size = 0.40, + show_polygon = TRUE, + use_overlap = FALSE, + polygon_feat_type = 'hex100', + polygon_bg_color = NA, + polygon_color = 'white', + polygon_line_size = 0.05, + jitter = c(25,25)) + +``` + +```{r, echo=FALSE, out.width="100%", fig.align="center", fig.cap="xxx"} +knitr::include_graphics("img/02_session2/23-spatInSituPlotPoints.png") +``` + + ```{r, eval=FALSE} hexbin25 <- tessellate(extent = ext(visiumHD_subset@feat_info$rna), shape = 'hexagon', @@ -751,17 +798,44 @@ visiumHD_subset = setPolygonInfo(gobject = visiumHD_subset, name = 'hex25', initialize = T) +showGiottoSpatialInfo(visiumHD_subset) + visiumHD_subset = addSpatialCentroidLocations(gobject = visiumHD_subset, poly_info = 'hex25') activeSpatUnit(visiumHD_subset) <- 'hex25' + +spatInSituPlotPoints(visiumHD_subset, + show_image = F, + feats = list('rna' = selected_feats), + feats_color_code = my_colors, + show_legend = F, + spat_unit = 'hex25', + point_size = 0.40, + show_polygon = TRUE, + use_overlap = FALSE, + polygon_feat_type = 'hex25', + polygon_bg_color = NA, + polygon_color = 'white', + polygon_line_size = 0.05, + jitter = c(25,25)) +``` + +```{r, echo=FALSE, out.width="100%", fig.align="center", fig.cap="xxx"} +knitr::include_graphics("img/02_session2/24-spatInSituPlotPoints.png") ``` + + + ```{r, eval=FALSE} visiumHD_subset = calculateOverlap(visiumHD_subset, spatial_info = 'hex25', feat_info = 'rna') +showGiottoSpatialInfo(visiumHD_subset) + + # convert overlap results to bin by gene matrix visiumHD_subset = overlapToMatrix(visiumHD_subset, poly_info = 'hex25', @@ -770,14 +844,234 @@ visiumHD_subset = overlapToMatrix(visiumHD_subset, visiumHD_subset <- filterGiotto(visiumHD_subset, expression_threshold = 1, - feat_det_in_min_cells = 25, - min_det_feats_per_cell = 25) + feat_det_in_min_cells = 3, + min_det_feats_per_cell = 5) + +activeSpatUnit(visiumHD_subset) + # normalize visiumHD_subset <- normalizeGiotto(visiumHD_subset, scalefactor = 1000, verbose = T) # add statistics visiumHD_subset <- addStatistics(visiumHD_subset) + +feature_data = fDataDT(visiumHD_subset) + +visiumHD_subset <- calculateHVF(visiumHD_subset, zscore_threshold = 1) + +``` + +### projections + +#### dimension reduction with projections +```{r, eval=FALSE} + +n_25_percent <- round(length(spatIDs(visiumHD_subset, 'hex25')) * 0.25) + +# pca projection on subset +visiumHD_subset <- runPCAprojection( + gobject = visiumHD_subset, + spat_unit = "hex25", + feats_to_use = 'hvf', + name = 'pca.projection', + set_seed = TRUE, + seed_number = 12345, + random_subset = n_25_percent +) + +showGiottoDimRed(visiumHD_subset) +plotPCA(visiumHD_subset, dim_reduction_name = 'pca.projection') +``` + +```{r, echo=FALSE, out.width="100%", fig.align="center", fig.cap="xxx"} +knitr::include_graphics("img/02_session2/25-PCA.png") +``` + + +```{r, eval=FALSE} +# umap projection on subset +visiumHD_subset <- runUMAPprojection( + gobject = visiumHD_subset, + spat_unit = "hex25", + dim_reduction_to_use = 'pca', + dim_reduction_name = "pca.projection", + dimensions_to_use = 1:10, + name = "umap.projection", + random_subset = n_25_percent, + n_neighbors = 10, + min_dist = 0.005, + n_threads = 4 +) + +showGiottoDimRed(visiumHD_subset) + +# plot UMAP, coloring cells/points based on nr_feats +plotUMAP(gobject = visiumHD_subset, + point_size = 1, + dim_reduction_name = 'umap.projection') + +``` + +```{r, echo=FALSE, out.width="100%", fig.align="center", fig.cap="xxx"} +knitr::include_graphics("img/02_session2/26-UMAP.png") +``` + + + +#### clustering with projection + +```{r, eval=FALSE} + +# subset to smaller giotto object +set.seed(1234) +subset_IDs = sample(x = spatIDs(visiumHD_subset, 'hex25'), size = n_25_percent) + +temp_gobject = subsetGiotto( + gobject = visiumHD_subset, + spat_unit = 'hex25', + cell_ids = subset_IDs +) + + +# hierarchical clustering +temp_gobject = doHclust(gobject = temp_gobject, + spat_unit = 'hex25', + k = 8, name = 'sub_hclust', + dim_reduction_to_use = 'pca', + dim_reduction_name = 'pca.projection', + dimensions_to_use = 1:10) + + +# show umap +dimPlot2D( + gobject = temp_gobject, + point_size = 2.5, + spat_unit = 'hex25', + dim_reduction_to_use = 'umap', + dim_reduction_name = 'umap.projection', + cell_color = 'sub_hclust' +) + +``` + +```{r, echo=FALSE, out.width="100%", fig.align="center", fig.cap="xxx"} +knitr::include_graphics("img/02_session2/27-dimPlot2D.png") +``` + + + +```{r, eval=FALSE} +# project clusterings back to full dataset +visiumHD_subset <- doClusterProjection( + target_gobject = visiumHD_subset, + source_gobject = temp_gobject, + spat_unit = "hex25", + source_cluster_labels = "sub_hclust", + reduction_method = 'pca', + reduction_name = 'pca.projection', + prob = FALSE, + knn_k = 5, + dimensions_to_use = 1:10 +) + +pDataDT(visiumHD_subset) + +dimPlot2D( + gobject = visiumHD_subset, + point_size = 1.5, + spat_unit = 'hex25', + dim_reduction_to_use = 'umap', + dim_reduction_name = 'umap.projection', + cell_color = 'knn_labels' +) + +``` + +```{r, echo=FALSE, out.width="100%", fig.align="center", fig.cap="xxx"} +knitr::include_graphics("img/02_session2/28-dimPlot2D.png") +``` + + + +```{r, eval=FALSE} +spatInSituPlotPoints(visiumHD_subset, + show_image = F, + feats = NULL, + show_legend = F, + spat_unit = 'hex25', + point_size = 0.5, + show_polygon = TRUE, + use_overlap = FALSE, + polygon_feat_type = 'hex25', + polygon_fill_as_factor = TRUE, + polygon_fill = 'knn_labels', + polygon_color = 'black', + polygon_line_size = 0.3) +``` + +```{r, echo=FALSE, out.width="100%", fig.align="center", fig.cap="xxx"} +knitr::include_graphics("img/02_session2/29-spatInSituPlotPoints.png") ``` +### niche clusters + + +```{r, eval=FALSE} +visiumHD_subset = createSpatialNetwork(visiumHD_subset, + name = 'kNN_network', + spat_unit = 'hex25', + method = 'kNN', + k = 6) + +pDataDT(visiumHD_subset) +visiumHD_subset = calculateSpatCellMetadataProportions(gobject = visiumHD_subset, + spat_unit = 'hex25', + feat_type = 'rna', + metadata_column = 'knn_labels', + spat_network = 'kNN_network') + +prop_table = getSpatialEnrichment(visiumHD_subset, name = 'proportion', output = 'data.table') +prop_matrix = GiottoUtils:::dt_to_matrix(prop_table) + +set.seed(1234) +prop_kmeans = kmeans(x = prop_matrix, centers = 10, iter.max = 1000, nstart = 100) +prop_kmeansDT = data.table::data.table(cell_ID = names(prop_kmeans$cluster), niche = prop_kmeans$cluster) + +visiumHD_subset = addCellMetadata(visiumHD_subset, + new_metadata = prop_kmeansDT, + by_column = T, + column_cell_ID = 'cell_ID') +pDataDT(visiumHD_subset) + +spatInSituPlotPoints(visiumHD_subset, + show_image = F, + feats = NULL, + show_legend = F, + spat_unit = 'hex25', + point_size = 0.5, + show_polygon = TRUE, + use_overlap = FALSE, + polygon_feat_type = 'hex25', + polygon_fill_as_factor = TRUE, + polygon_fill = 'niche', + polygon_color = 'black', + polygon_line_size = 0.3) + +``` + +```{r, echo=FALSE, out.width="100%", fig.align="center", fig.cap="xxx"} +knitr::include_graphics("img/02_session2/30-spatInSituPlotPoints.png") +``` + + + ## Database backend - Work in progress, but coming soon! + +Memory problems: +- data ingestion +- spatial operations +- matrix operations +- matrix and spatial geometry object sizes + + diff --git a/docs/404.html b/docs/404.html index b85d1e4..644c3ee 100644 --- a/docs/404.html +++ b/docs/404.html @@ -325,6 +325,8 @@
  • 9.6 Hexbin 25
  • 9.7 Database backend - Work in progress, but coming soon!
    • diff --git a/docs/img/02_session2/0-spatInSituPlotPoints.png b/docs/img/02_session2/0-spatInSituPlotPoints.png index 83a6f19..4832d47 100644 Binary files a/docs/img/02_session2/0-spatInSituPlotPoints.png and b/docs/img/02_session2/0-spatInSituPlotPoints.png differ diff --git a/docs/img/02_session2/1-spatInSituPlotPoints.png b/docs/img/02_session2/1-spatInSituPlotPoints.png index 58ddb68..a769766 100644 Binary files a/docs/img/02_session2/1-spatInSituPlotPoints.png and b/docs/img/02_session2/1-spatInSituPlotPoints.png differ diff --git a/docs/img/02_session2/12-UMAP.png b/docs/img/02_session2/12-UMAP.png index 3e6e878..bbbb2ae 100644 Binary files a/docs/img/02_session2/12-UMAP.png and b/docs/img/02_session2/12-UMAP.png differ diff --git a/docs/img/02_session2/21-spatInSituPlotPoints.png b/docs/img/02_session2/21-spatInSituPlotPoints.png index 47dcfb8..96f4ac4 100644 Binary files a/docs/img/02_session2/21-spatInSituPlotPoints.png and b/docs/img/02_session2/21-spatInSituPlotPoints.png differ diff --git a/docs/img/02_session2/7-UMAP.png b/docs/img/02_session2/7-UMAP.png index 2c2cd83..95215d1 100644 Binary files a/docs/img/02_session2/7-UMAP.png and b/docs/img/02_session2/7-UMAP.png differ diff --git a/docs/img/02_session2/9-spatInSituPlotPoints.png b/docs/img/02_session2/9-spatInSituPlotPoints.png index cb11008..723b9d1 100644 Binary files a/docs/img/02_session2/9-spatInSituPlotPoints.png and b/docs/img/02_session2/9-spatInSituPlotPoints.png differ diff --git a/docs/reference-keys.txt b/docs/reference-keys.txt index b792810..86c8074 100644 --- a/docs/reference-keys.txt +++ b/docs/reference-keys.txt @@ -506,3 +506,15 @@ fig:unnamed-chunk-61 fig:unnamed-chunk-63 fig:unnamed-chunk-65 fig:unnamed-chunk-67 +fig:unnamed-chunk-74 +fig:unnamed-chunk-76 +fig:unnamed-chunk-79 +fig:unnamed-chunk-81 +fig:unnamed-chunk-83 +fig:unnamed-chunk-85 +fig:unnamed-chunk-87 +fig:unnamed-chunk-89 +projections +dimension-reduction-with-projections +clustering-with-projection +niche-clusters diff --git a/docs/search_index.json b/docs/search_index.json index 09f7d89..a766d44 100644 --- a/docs/search_index.json +++ b/docs/search_index.json @@ -1 +1 @@ -[["visium-hd.html", "9 Visium HD 9.1 Objective 9.2 Background 9.3 Data Ingestion 9.4 Giotto object at Hexbin 400 9.5 Hexbin 100 9.6 Hexbin 25 9.7 Database backend - Work in progress, but coming soon!", " 9 Visium HD Ruben Dries & Edward C. Ruiz August 6th 2024 9.1 Objective This tutorial demonstrates how to process Visium HD data at the highest 2 micron bin resolution using Giotto Suite. Notably, a similar strategy can be used for other spatial sequencing methods that operate at the subcellular level, including: - Stereo-seq - Seq-Scope - Open-ST The resulting datasets from all these technologies can be very large since they provide both a high spatial resolution and genome-wide capture of all transcripts. 9.2 Background 9.2.1 Visium HD Technology Figure 9.1: Overview of Visium HD. Source: 10X Genomics Visium HD is a spatial transcriptomics technology recently developed by 10X Genomics. Details about this platform are discussed on the official 10X Genomics Visium HD website and the preprint by Oliveira et al. 2024 on bioRxiv. Visium HD has a 2 micron bin size resolution. The default SpaceRanger pipeline from 10X Genomics also returns aggregated data at the 8 and 16 micron bin size. 9.2.2 Colorectal Cancer Sample Figure 9.2: Colorectal Cancer Overview. Source: 10X Genomics For this tutorial we will be using the publicly available Colorectal Cancer Visium HD dataset. Details about this dataset and a link to download the raw data can be found at the 10X Genomics website. 9.3 Data Ingestion 9.3.1 Visium HD output data format Figure 9.3: File structure of Visium HD data processed with spaceranger pipeline. Visium HD data processed with the spaceranger pipeline is organized in this format containing various files associated with the sample. The files highlighted in yellow are what we will be using to read in these datasets. 9.3.2 Mini Visium HD dataset library(Giotto) # set up paths data_path <- "data/02_session2/" save_dir <- "results/02_session2/" dir.create(save_dir, recursive = TRUE) # download the mini dataset and untar options("timeout" = Inf) download.file( url = "https://zenodo.org/records/13212855/files/workshop_VisiumHD.zip?download=1", destfile = file.path(save_dir, "workshop_visiumHD.zip") ) untar(tarfile = file.path(save_dir, "workshop_visiumHD.zip"), exdir = data_path) 9.3.3 Giotto Visium HD convenience function The easiest way to read in Visium HD data in Giotto is through our convenience function. This function will automatically read in the data at your desired resolution, align the images, and finally create a Giotto Object. # importVisiumHD() 9.3.4 Read in data manually However, for this tutorial we will illustrate how to create your own Giotto object in a step-by-step manner, which can also be applied to other similar technologies as discussed in the Objective section. 9.3.4.1 Raw expression data expression_path <- file.path(data_path, 'workshop_visiumHD/Human_Colorectal_Cancer_workshop/square_002um/raw_feature_bc_matrix') expr_results <- get10Xmatrix(path_to_data = expression_path, gene_column_index = 1) 9.3.4.2 Tissue positions data tissue_positions_path <- file.path(data_path, 'workshop_visiumHD/Human_Colorectal_Cancer_workshop/square_002um/spatial/tissue_positions.parquet') tissue_positions <- data.table::as.data.table(arrow::read_parquet(tissue_positions_path)) 9.3.4.3 Merge expression and 2 micron position data # convert expression matrix to minimal data.frame or data.table object matrix_tile_dt <- data.table::as.data.table(Matrix::summary(expr_results)) genes <- expr_results@Dimnames[[1]] samples <- expr_results@Dimnames[[2]] matrix_tile_dt[, gene := genes[i]] matrix_tile_dt[, pixel := samples[j]] Figure 9.4: Genes expressed for each 2 µm pixel. # merge data.table matrix and spatial coordinates to create input for Giotto Polygons expr_pos_data <- data.table::merge.data.table(matrix_tile_dt, tissue_positions, by.x = 'pixel', by.y = 'barcode') expr_pos_data <- expr_pos_data[,.(pixel, pxl_row_in_fullres, pxl_col_in_fullres, gene, x)] colnames(expr_pos_data) = c('pixel', 'x', 'y', 'gene', 'count') Figure 9.5: Genes expressed with count for each 2 µm pixel. 9.4 Giotto object at Hexbin 400 9.4.1 create giotto points giotto_points = createGiottoPoints(x = expr_pos_data[,.(x, y, gene, pixel, count)]) 9.4.2 create giotto polygons 9.4.2.1 Tiling and aggregation The Visium HD data is organized in a grid format. We can aggregate the data into larger bins to reduce the dimensionality of the data. Giotto Suite provides options to bin data not only with squares, but also through hexagons and triangles. Here we use a hexagon tesselation to aggregate the data into arbitrary cells. # create giotto polygons, here we create hexagons hexbin400 <- tessellate(extent = ext(giotto_points), shape = 'hexagon', shape_size = 400, name = 'hex400') plot(hexbin400) Figure 9.6: Giotto polygon in a hexagon shape for overlapping visium HD expression data. 9.4.3 combine Giotto points and polygons to create Giotto object instrs = createGiottoInstructions( save_dir = save_dir, save_plot = TRUE, show_plot = FALSE, return_plot = FALSE ) visiumHD = createGiottoObjectSubcellular(gpoints = list('rna' = giotto_points), gpolygons = list('hex400' = hexbin400), instructions = instrs) visiumHD = addSpatialCentroidLocations(gobject = visiumHD, poly_info = 'hex400') Visualize the Giotto object. Make sure to set expand_counts = TRUE to expand the counts column. Each spatial bin can have multiple transcripts/UMIs. This is different compared to in situ technologies like seqFISH, MERFISH, Nanostring CosMx or Xenium. feature_data = fDataDT(visiumHD) spatInSituPlotPoints(visiumHD, show_image = F, feats = list('rna' = feature_data$feat_ID[10:20]), show_legend = T, spat_unit = 'hex400', point_size = 0.25, show_polygon = TRUE, use_overlap = FALSE, polygon_feat_type = 'hex400', polygon_bg_color = NA, polygon_color = 'white', polygon_line_size = 0.1, expand_counts = TRUE, count_info_column = 'count', jitter = c(25,25)) Figure 9.7: Overlap of gene expression with the hex400 polygons. Each dot represents a single gene. Jitter used to better vizualize individual genes. You can set plot_method = scattermore or scattermost to convert high-resolution images to low(er) resolution rasterized images. It’s usually faster and will save on disk space. spatInSituPlotPoints(visiumHD, show_image = F, feats = list('rna' = feature_data$feat_ID[10:20]), show_legend = T, spat_unit = 'hex400', point_size = 0.25, show_polygon = TRUE, use_overlap = FALSE, polygon_feat_type = 'hex400', polygon_bg_color = NA, polygon_color = 'white', polygon_line_size = 0.1, expand_counts = TRUE, count_info_column = 'count', jitter = c(25,25), plot_method = 'scattermore') 9.4.4 Process Giotto object 9.4.4.1 calculate overlap between points and polygons # calculate overlap between points and polygons visiumHD = calculateOverlap(visiumHD, spatial_info = 'hex400', feat_info = 'rna') showGiottoSpatialInfo(visiumHD) 9.4.4.2 convert overlap results to a bin-by-gene matrix # convert overlap results to bin by gene matrix visiumHD = overlapToMatrix(visiumHD, poly_info = 'hex400', feat_info = 'rna', name = 'raw') activeSpatUnit(visiumHD) 9.4.4.3 default processing steps visiumHD <- filterGiotto(visiumHD, expression_threshold = 1, feat_det_in_min_cells = 5, min_det_feats_per_cell = 25) # normalize visiumHD <- normalizeGiotto(visiumHD, scalefactor = 1000, verbose = T) # add statistics visiumHD <- addStatistics(visiumHD) 9.4.4.4 visualize number of features At the centroid level. # each dot here represents a 200x200 aggregation of spatial barcodes (bin size 200) spatPlot2D(gobject = visiumHD, cell_color = "nr_feats", color_as_factor = F, point_size = 2.5) Figure 9.8: Number of features detected in each of the centroids. Using the spatial polygon (hexagon) tiles spatInSituPlotPoints(visiumHD, show_image = F, feats = NULL, show_legend = F, spat_unit = 'hex400', point_size = 0.1, show_polygon = TRUE, use_overlap = TRUE, polygon_feat_type = 'hex400', polygon_fill = 'nr_feats', polygon_fill_as_factor = F, polygon_bg_color = NA, polygon_color = 'white', polygon_line_size = 0.1) Figure 9.9: Number of features detected in each of the hex400 polygons. 9.4.5 Dimension reduction + clustering 9.4.5.1 Highly variable features + PCA visiumHD <- calculateHVF(visiumHD, zscore_threshold = 1) visiumHD <- runPCA(visiumHD, expression_values = 'normalized', feats_to_use = 'hvf') screePlot(visiumHD, ncp = 30) plotPCA(visiumHD) 9.4.5.2 UMAP reduction for visualization visiumHD <- runUMAP(visiumHD, dimensions_to_use = 1:14, n_threads = 10) plotUMAP(gobject = visiumHD, point_size = 1) 9.4.5.3 Create network based on expression similarity + graph partition cluster # sNN network (default) visiumHD <- createNearestNetwork(visiumHD, dimensions_to_use = 1:14, k = 5) ## leiden clustering #### visiumHD <- doLeidenClusterIgraph(visiumHD, resolution = 0.5, n_iterations = 1000, spat_unit = 'hex400') plotUMAP(gobject = visiumHD, cell_color = 'leiden_clus', point_size = 1.5, show_NN_network = F, edge_alpha = 0.05) Figure 9.10: Leiden clustering for the hex400 bins. spatInSituPlotPoints(visiumHD, show_image = F, feats = NULL, show_legend = F, spat_unit = 'hex400', point_size = 0.25, show_polygon = TRUE, use_overlap = FALSE, polygon_feat_type = 'hex400', polygon_fill_as_factor = TRUE, polygon_fill = 'leiden_clus', polygon_color = 'black', polygon_line_size = 0.3) Figure 9.11: Spat plot for hex400 bin colored by leiden clusters. 9.5 Hexbin 100 Goal is to create a higher resolution bin (hex100) and add to the Giotto object. 9.5.1 Standard subcellular pipeline hexbin100 <- tessellate(extent = ext(visiumHD), shape = 'hexagon', shape_size = 100, name = 'hex100') visiumHD = setPolygonInfo(gobject = visiumHD, x = hexbin100, name = 'hex100', initialize = T) visiumHD = addSpatialCentroidLocations(gobject = visiumHD, poly_info = 'hex100') Set active spatial unit. This can also be set manually. activeSpatUnit(visiumHD) <- 'hex100' spatInSituPlotPoints(visiumHD, show_image = F, feats = list('rna' = feature_data$feat_ID[1:20]), show_legend = T, spat_unit = 'hex100', point_size = 0.1, show_polygon = TRUE, use_overlap = FALSE, polygon_feat_type = 'hex100', polygon_bg_color = NA, polygon_color = 'white', polygon_line_size = 0.2, expand_counts = TRUE, count_info_column = 'count', jitter = c(25,25)) Figure 9.12: Polygon overlay of hex100 bins over 2 µm pixel. Jitter applied to vizualize individual features. visiumHD = calculateOverlap(visiumHD, spatial_info = 'hex100', feat_info = 'rna') visiumHD = overlapToMatrix(visiumHD, poly_info = 'hex100', feat_info = 'rna', name = 'raw') visiumHD <- filterGiotto(visiumHD, expression_threshold = 1, feat_det_in_min_cells = 10, min_det_feats_per_cell = 10) visiumHD <- normalizeGiotto(visiumHD, scalefactor = 1000, verbose = T) visiumHD <- addStatistics(visiumHD) pDataDT(visiumHD, spat_unit = 'hex100') pDataDT(visiumHD, spat_unit = 'hex400') ## dimension reduction #### # --------------------------- # visiumHD <- calculateHVF(visiumHD, zscore_threshold = 1) visiumHD <- runPCA(visiumHD, expression_values = 'normalized', feats_to_use = 'hvf') plotPCA(visiumHD) visiumHD <- runUMAP(visiumHD, dimensions_to_use = 1:14, n_threads = 10) # plot UMAP, coloring cells/points based on nr_feats plotUMAP(gobject = visiumHD, point_size = 2) Figure 9.13: UMAP for the hex100 bin. # sNN network (default) visiumHD <- createNearestNetwork(visiumHD, dimensions_to_use = 1:14, k = 5) ## leiden clustering #### visiumHD <- doLeidenClusterIgraph(visiumHD, resolution = 0.2, n_iterations = 1000) plotUMAP(gobject = visiumHD, cell_color = 'leiden_clus', point_size = 1.5, show_NN_network = F, edge_alpha = 0.05) Figure 9.14: UMAP for the hex100 bin colored by ledien clusters. spatInSituPlotPoints(visiumHD, show_image = F, feats = NULL, show_legend = F, spat_unit = 'hex100', point_size = 0.5, show_polygon = TRUE, use_overlap = FALSE, polygon_feat_type = 'hex100', polygon_fill_as_factor = TRUE, polygon_fill = 'leiden_clus', polygon_color = 'black', polygon_line_size = 0.3) Figure 9.15: Spat plot for the hex100 bin colored by leiden clusters. 9.5.2 Spatial expression patterns 9.5.2.1 Identify single genes featData = fDataDT(visiumHD) hvf_genes = featData[hvf == 'yes']$feat_ID visiumHD = createSpatialNetwork(visiumHD, name = 'kNN_network', spat_unit = 'hex100', method = 'kNN', k = 8) ranktest = binSpect(visiumHD, spat_unit = 'hex100', subset_feats = hvf_genes, bin_method = 'rank', calc_hub = FALSE, do_fisher_test = TRUE, spatial_network_name = 'kNN_network') Visualize top 2 ranked spatial genes per expression bin: set0 = ranktest[high_expr < 50][1:2]$feats set1 = ranktest[high_expr > 50 & high_expr < 100][1:2]$feats set2 = ranktest[high_expr > 100 & high_expr < 200][1:2]$feats set3 = ranktest[high_expr > 200 & high_expr < 400][1:2]$feats set4 = ranktest[high_expr > 400 & high_expr < 1000][1:2]$feats set5 = ranktest[high_expr > 1000][1:2]$feats spatFeatPlot2D(visiumHD, expression_values = 'scaled', feats = c(set0, set1, set2), gradient_style = "sequential", cell_color_gradient = c('blue', 'white', 'yellow', 'orange', 'red', 'darkred'), cow_n_col = 2, point_size = 1) Figure 9.16: Spat feature plot showing gene expression for the top 2 ranked spatial genes per expression bin (<50, >50 and >100) across the hex100 bin. spatFeatPlot2D(visiumHD, expression_values = 'scaled', feats = c(set3, set4, set5), gradient_style = "sequential", cell_color_gradient = c('blue', 'white', 'yellow', 'orange', 'red', 'darkred'), cow_n_col = 2, point_size = 1) Figure 9.17: Spat feature plot showing gene expression for the top 2 ranked spatial genes per expression bin (>200, >400 and >1000) across the hex100 bin. 9.5.2.2 Spatial co-expression modules ext_spatial_genes = ranktest[adj.p.value < 0.001]$feats spat_cor_netw_DT = detectSpatialCorFeats(visiumHD, method = 'network', spatial_network_name = 'kNN_network', subset_feats = ext_spatial_genes) # cluster spatial genes spat_cor_netw_DT = clusterSpatialCorFeats(spat_cor_netw_DT, name = 'spat_netw_clus', k = 16) # visualize clusters heatmSpatialCorFeats(visiumHD, spatCorObject = spat_cor_netw_DT, use_clus_name = 'spat_netw_clus', heatmap_legend_param = list(title = NULL)) Figure 9.18: Heatmap showing spatially correlated genes split into 16 clusters. # create metagene enrichment score for clusters cluster_genes_DT = showSpatialCorFeats(spat_cor_netw_DT, use_clus_name = 'spat_netw_clus', show_top_feats = 1) cluster_genes = cluster_genes_DT$clus; names(cluster_genes) = cluster_genes_DT$feat_ID visiumHD = createMetafeats(visiumHD, expression_values = 'normalized', feat_clusters = cluster_genes, name = 'cluster_metagene') showGiottoSpatEnrichments(visiumHD) spatCellPlot(visiumHD, spat_enr_names = 'cluster_metagene', gradient_style = "sequential", cell_color_gradient = c('blue', 'white', 'yellow', 'orange', 'red', 'darkred'), cell_annotation_values = as.character(c(1:4)), point_size = 1, cow_n_col = 2) Figure 9.19: Spat plot vizualizing metagenes (1-4) based on spatially correlated genes vizualized on the hex100 bin spatCellPlot(visiumHD, spat_enr_names = 'cluster_metagene', gradient_style = "sequential", cell_color_gradient = c('blue', 'white', 'yellow', 'orange', 'red', 'darkred'), cell_annotation_values = as.character(c(5:8)), point_size = 1, cow_n_col = 2) Figure 9.20: Spat plot vizualizing metagenes (5-8) based on spatially correlated genes vizualized on the hex100 bin spatCellPlot(visiumHD, spat_enr_names = 'cluster_metagene', gradient_style = "sequential", cell_color_gradient = c('blue', 'white', 'yellow', 'orange', 'red', 'darkred'), cell_annotation_values = as.character(c(9:12)), point_size = 1, cow_n_col = 2) Figure 9.21: Spat plot vizualizing metagenes (9-12) based on spatially correlated genes vizualized on the hex100 bin spatCellPlot(visiumHD, spat_enr_names = 'cluster_metagene', gradient_style = "sequential", cell_color_gradient = c('blue', 'white', 'yellow', 'orange', 'red', 'darkred'), cell_annotation_values = as.character(c(13:16)), point_size = 1, cow_n_col = 2) Figure 9.22: Spat plot vizualizing metagenes (13-16) based on spatially correlated genes vizualized on the hex100 bin 9.5.2.3 Plot spatial gene groups balanced_genes = getBalancedSpatCoexpressionFeats(spatCorObject = spat_cor_netw_DT, maximum = 5) selected_feats = names(balanced_genes) # give genes from same cluster same color distinct_colors = getDistinctColors(n = 20) names(distinct_colors) = 1:20 my_colors = distinct_colors[balanced_genes] names(my_colors) = names(balanced_genes) spatInSituPlotPoints(visiumHD, show_image = F, feats = list('rna' = selected_feats), feats_color_code = my_colors, show_legend = F, spat_unit = 'hex100', point_size = 0.20, show_polygon = FALSE, use_overlap = FALSE, polygon_feat_type = 'hex100', polygon_bg_color = NA, polygon_color = 'white', polygon_line_size = 0.01, expand_counts = TRUE, count_info_column = 'count', jitter = c(25,25)) Figure 9.23: Coloring individual features based on the spatially correlated gene clusters. 9.6 Hexbin 25 Goal is to create a higher resolution bin (hex25) and add to the Giotto object. We will aim to identify individual cell types and local neighborhood niches. 9.6.1 Subcellular workflow and projection functions filter and normalization workflow PCA projection visiumHD_subset = subsetGiottoLocs(gobject = visiumHD, x_min = 16000, x_max = 20000, y_min = 44250, y_max = 45500) hexbin25 <- tessellate(extent = ext(visiumHD_subset@feat_info$rna), shape = 'hexagon', shape_size = 25, name = 'hex25') visiumHD_subset = setPolygonInfo(gobject = visiumHD_subset, x = hexbin25, name = 'hex25', initialize = T) visiumHD_subset = addSpatialCentroidLocations(gobject = visiumHD_subset, poly_info = 'hex25') activeSpatUnit(visiumHD_subset) <- 'hex25' visiumHD_subset = calculateOverlap(visiumHD_subset, spatial_info = 'hex25', feat_info = 'rna') # convert overlap results to bin by gene matrix visiumHD_subset = overlapToMatrix(visiumHD_subset, poly_info = 'hex25', feat_info = 'rna', name = 'raw') visiumHD_subset <- filterGiotto(visiumHD_subset, expression_threshold = 1, feat_det_in_min_cells = 25, min_det_feats_per_cell = 25) # normalize visiumHD_subset <- normalizeGiotto(visiumHD_subset, scalefactor = 1000, verbose = T) # add statistics visiumHD_subset <- addStatistics(visiumHD_subset) 9.7 Database backend - Work in progress, but coming soon! "],["404.html", "Page not found", " Page not found The page you requested cannot be found (perhaps it was moved or renamed). You may want to try searching to find the page's new location, or use the table of contents to find the page you are looking for. "]] +[["visium-hd.html", "9 Visium HD 9.1 Objective 9.2 Background 9.3 Data Ingestion 9.4 Giotto object at Hexbin 400 9.5 Hexbin 100 9.6 Hexbin 25 9.7 Database backend - Work in progress, but coming soon!", " 9 Visium HD Ruben Dries & Edward C. Ruiz August 6th 2024 9.1 Objective This tutorial demonstrates how to process Visium HD data at the highest 2 micron bin resolution using Giotto Suite. Notably, a similar strategy can be used for other spatial sequencing methods that operate at the subcellular level, including: - Stereo-seq - Seq-Scope - Open-ST The resulting datasets from all these technologies can be very large since they provide both a high spatial resolution and genome-wide capture of all transcripts. 9.2 Background 9.2.1 Visium HD Technology Figure 9.1: Overview of Visium HD. Source: 10X Genomics Visium HD is a spatial transcriptomics technology recently developed by 10X Genomics. Details about this platform are discussed on the official 10X Genomics Visium HD website and the preprint by Oliveira et al. 2024 on bioRxiv. Visium HD has a 2 micron bin size resolution. The default SpaceRanger pipeline from 10X Genomics also returns aggregated data at the 8 and 16 micron bin size. 9.2.2 Colorectal Cancer Sample Figure 9.2: Colorectal Cancer Overview. Source: 10X Genomics For this tutorial we will be using the publicly available Colorectal Cancer Visium HD dataset. Details about this dataset and a link to download the raw data can be found at the 10X Genomics website. 9.3 Data Ingestion 9.3.1 Visium HD output data format Figure 9.3: File structure of Visium HD data processed with spaceranger pipeline. Visium HD data processed with the spaceranger pipeline is organized in this format containing various files associated with the sample. The files highlighted in yellow are what we will be using to read in these datasets. 9.3.2 Mini Visium HD dataset library(Giotto) # set up paths data_path <- "data/02_session2/" save_dir <- "results/02_session2/" dir.create(save_dir, recursive = TRUE) # download the mini dataset and untar options("timeout" = Inf) download.file( url = "https://zenodo.org/records/13226158/files/workshop_VisiumHD.zip?download=1", destfile = file.path(save_dir, "workshop_visiumHD.zip") ) untar(tarfile = file.path(save_dir, "workshop_visiumHD.zip"), exdir = data_path) 9.3.3 Giotto Visium HD convenience function The easiest way to read in Visium HD data in Giotto is through our convenience function. This function will automatically read in the data at your desired resolution, align the images, and finally create a Giotto Object. # importVisiumHD() 9.3.4 Read in data manually However, for this tutorial we will illustrate how to create your own Giotto object in a step-by-step manner, which can also be applied to other similar technologies as discussed in the Objective section. 9.3.4.1 Raw expression data expression_path <- file.path(data_path, '/Human_Colorectal_Cancer_workshop/square_002um/raw_feature_bc_matrix') expr_results <- get10Xmatrix(path_to_data = expression_path, gene_column_index = 1) 9.3.4.2 Tissue positions data tissue_positions_path <- file.path(data_path, '/Human_Colorectal_Cancer_workshop/square_002um/spatial/tissue_positions.parquet') tissue_positions <- data.table::as.data.table(arrow::read_parquet(tissue_positions_path)) 9.3.4.3 Merge expression and 2 micron position data # convert expression matrix to minimal data.frame or data.table object matrix_tile_dt <- data.table::as.data.table(Matrix::summary(expr_results)) genes <- expr_results@Dimnames[[1]] samples <- expr_results@Dimnames[[2]] matrix_tile_dt[, gene := genes[i]] matrix_tile_dt[, pixel := samples[j]] Figure 9.4: Genes expressed for each 2 µm pixel. # merge data.table matrix and spatial coordinates to create input for Giotto Polygons expr_pos_data <- data.table::merge.data.table(matrix_tile_dt, tissue_positions, by.x = 'pixel', by.y = 'barcode') expr_pos_data <- expr_pos_data[,.(pixel, pxl_row_in_fullres, pxl_col_in_fullres, gene, x)] colnames(expr_pos_data) = c('pixel', 'x', 'y', 'gene', 'count') Figure 9.5: Genes expressed with count for each 2 µm pixel. 9.4 Giotto object at Hexbin 400 9.4.1 create giotto points giotto_points = createGiottoPoints(x = expr_pos_data[,.(x, y, gene, pixel, count)]) 9.4.2 create giotto polygons 9.4.2.1 Tiling and aggregation The Visium HD data is organized in a grid format. We can aggregate the data into larger bins to reduce the dimensionality of the data. Giotto Suite provides options to bin data not only with squares, but also through hexagons and triangles. Here we use a hexagon tesselation to aggregate the data into arbitrary cells. # create giotto polygons, here we create hexagons hexbin400 <- tessellate(extent = ext(giotto_points), shape = 'hexagon', shape_size = 400, name = 'hex400') plot(hexbin400) Figure 9.6: Giotto polygon in a hexagon shape for overlapping visium HD expression data. 9.4.3 combine Giotto points and polygons to create Giotto object instrs = createGiottoInstructions( save_dir = save_dir, save_plot = TRUE, show_plot = FALSE, return_plot = FALSE ) visiumHD = createGiottoObjectSubcellular(gpoints = list('rna' = giotto_points), gpolygons = list('hex400' = hexbin400), instructions = instrs) visiumHD = addSpatialCentroidLocations(gobject = visiumHD, poly_info = 'hex400') Visualize the Giotto object. Make sure to set expand_counts = TRUE to expand the counts column. Each spatial bin can have multiple transcripts/UMIs. This is different compared to in situ technologies like seqFISH, MERFISH, Nanostring CosMx or Xenium. feature_data = fDataDT(visiumHD) spatInSituPlotPoints(visiumHD, show_image = F, feats = list('rna' = feature_data$feat_ID[10:20]), show_legend = T, spat_unit = 'hex400', point_size = 0.25, show_polygon = TRUE, use_overlap = FALSE, polygon_feat_type = 'hex400', polygon_bg_color = NA, polygon_color = 'white', polygon_line_size = 0.1, expand_counts = TRUE, count_info_column = 'count', jitter = c(25,25)) Figure 9.7: Overlap of gene expression with the hex400 polygons. Each dot represents a single gene. Jitter used to better vizualize individual genes. You can set plot_method = scattermore or scattermost to convert high-resolution images to low(er) resolution rasterized images. It’s usually faster and will save on disk space. spatInSituPlotPoints(visiumHD, show_image = F, feats = list('rna' = feature_data$feat_ID[10:20]), show_legend = T, spat_unit = 'hex400', point_size = 0.25, show_polygon = TRUE, use_overlap = FALSE, polygon_feat_type = 'hex400', polygon_bg_color = NA, polygon_color = 'white', polygon_line_size = 0.1, expand_counts = TRUE, count_info_column = 'count', jitter = c(25,25), plot_method = 'scattermore') 9.4.4 Process Giotto object 9.4.4.1 calculate overlap between points and polygons # calculate overlap between points and polygons visiumHD = calculateOverlap(visiumHD, spatial_info = 'hex400', feat_info = 'rna') showGiottoSpatialInfo(visiumHD) 9.4.4.2 convert overlap results to a bin-by-gene matrix # convert overlap results to bin by gene matrix visiumHD = overlapToMatrix(visiumHD, poly_info = 'hex400', feat_info = 'rna', name = 'raw') activeSpatUnit(visiumHD) 9.4.4.3 default processing steps visiumHD <- filterGiotto(visiumHD, expression_threshold = 1, feat_det_in_min_cells = 5, min_det_feats_per_cell = 25) # normalize visiumHD <- normalizeGiotto(visiumHD, scalefactor = 1000, verbose = T) # add statistics visiumHD <- addStatistics(visiumHD) 9.4.4.4 visualize number of features At the centroid level. # each dot here represents a 200x200 aggregation of spatial barcodes (bin size 200) spatPlot2D(gobject = visiumHD, cell_color = "nr_feats", color_as_factor = F, point_size = 2.5) Figure 9.8: Number of features detected in each of the centroids. Using the spatial polygon (hexagon) tiles spatInSituPlotPoints(visiumHD, show_image = F, feats = NULL, show_legend = F, spat_unit = 'hex400', point_size = 0.1, show_polygon = TRUE, use_overlap = TRUE, polygon_feat_type = 'hex400', polygon_fill = 'nr_feats', polygon_fill_as_factor = F, polygon_bg_color = NA, polygon_color = 'white', polygon_line_size = 0.1) Figure 9.9: Number of features detected in each of the hex400 polygons. 9.4.5 Dimension reduction + clustering 9.4.5.1 Highly variable features + PCA visiumHD <- calculateHVF(visiumHD, zscore_threshold = 1) visiumHD <- runPCA(visiumHD, expression_values = 'normalized', feats_to_use = 'hvf') screePlot(visiumHD, ncp = 30) plotPCA(visiumHD) 9.4.5.2 UMAP reduction for visualization visiumHD <- runUMAP(visiumHD, dimensions_to_use = 1:14, n_threads = 10) plotUMAP(gobject = visiumHD, point_size = 1) 9.4.5.3 Create network based on expression similarity + graph partition cluster # sNN network (default) visiumHD <- createNearestNetwork(visiumHD, dimensions_to_use = 1:14, k = 5) ## leiden clustering #### visiumHD <- doLeidenClusterIgraph(visiumHD, resolution = 0.5, n_iterations = 1000, spat_unit = 'hex400') plotUMAP(gobject = visiumHD, cell_color = 'leiden_clus', point_size = 1.5, show_NN_network = F, edge_alpha = 0.05) Figure 9.10: Leiden clustering for the hex400 bins. spatInSituPlotPoints(visiumHD, show_image = F, feats = NULL, show_legend = F, spat_unit = 'hex400', point_size = 0.25, show_polygon = TRUE, use_overlap = FALSE, polygon_feat_type = 'hex400', polygon_fill_as_factor = TRUE, polygon_fill = 'leiden_clus', polygon_color = 'black', polygon_line_size = 0.3) Figure 9.11: Spat plot for hex400 bin colored by leiden clusters. 9.5 Hexbin 100 Goal is to create a higher resolution bin (hex100) and add to the Giotto object. 9.5.1 Standard subcellular pipeline hexbin100 <- tessellate(extent = ext(visiumHD), shape = 'hexagon', shape_size = 100, name = 'hex100') visiumHD = setPolygonInfo(gobject = visiumHD, x = hexbin100, name = 'hex100', initialize = T) visiumHD = addSpatialCentroidLocations(gobject = visiumHD, poly_info = 'hex100') Set active spatial unit. This can also be set manually. activeSpatUnit(visiumHD) <- 'hex100' spatInSituPlotPoints(visiumHD, show_image = F, feats = list('rna' = feature_data$feat_ID[1:20]), show_legend = T, spat_unit = 'hex100', point_size = 0.1, show_polygon = TRUE, use_overlap = FALSE, polygon_feat_type = 'hex100', polygon_bg_color = NA, polygon_color = 'white', polygon_line_size = 0.2, expand_counts = TRUE, count_info_column = 'count', jitter = c(25,25)) Figure 9.12: Polygon overlay of hex100 bins over 2 µm pixel. Jitter applied to vizualize individual features. visiumHD = calculateOverlap(visiumHD, spatial_info = 'hex100', feat_info = 'rna') visiumHD = overlapToMatrix(visiumHD, poly_info = 'hex100', feat_info = 'rna', name = 'raw') visiumHD <- filterGiotto(visiumHD, expression_threshold = 1, feat_det_in_min_cells = 10, min_det_feats_per_cell = 10) visiumHD <- normalizeGiotto(visiumHD, scalefactor = 1000, verbose = T) visiumHD <- addStatistics(visiumHD) pDataDT(visiumHD, spat_unit = 'hex100') pDataDT(visiumHD, spat_unit = 'hex400') ## dimension reduction #### # --------------------------- # visiumHD <- calculateHVF(visiumHD, zscore_threshold = 1) visiumHD <- runPCA(visiumHD, expression_values = 'normalized', feats_to_use = 'hvf') plotPCA(visiumHD) visiumHD <- runUMAP(visiumHD, dimensions_to_use = 1:14, n_threads = 10) # plot UMAP, coloring cells/points based on nr_feats plotUMAP(gobject = visiumHD, point_size = 2) Figure 9.13: UMAP for the hex100 bin. # sNN network (default) visiumHD <- createNearestNetwork(visiumHD, dimensions_to_use = 1:14, k = 5) ## leiden clustering #### visiumHD <- doLeidenClusterIgraph(visiumHD, resolution = 0.2, n_iterations = 1000) plotUMAP(gobject = visiumHD, cell_color = 'leiden_clus', point_size = 1.5, show_NN_network = F, edge_alpha = 0.05) Figure 9.14: UMAP for the hex100 bin colored by ledien clusters. spatInSituPlotPoints(visiumHD, show_image = F, feats = NULL, show_legend = F, spat_unit = 'hex100', point_size = 0.5, show_polygon = TRUE, use_overlap = FALSE, polygon_feat_type = 'hex100', polygon_fill_as_factor = TRUE, polygon_fill = 'leiden_clus', polygon_color = 'black', polygon_line_size = 0.3) Figure 9.15: Spat plot for the hex100 bin colored by leiden clusters. 9.5.2 Spatial expression patterns 9.5.2.1 Identify single genes featData = fDataDT(visiumHD) hvf_genes = featData[hvf == 'yes']$feat_ID visiumHD = createSpatialNetwork(visiumHD, name = 'kNN_network', spat_unit = 'hex100', method = 'kNN', k = 8) ranktest = binSpect(visiumHD, spat_unit = 'hex100', subset_feats = hvf_genes, bin_method = 'rank', calc_hub = FALSE, do_fisher_test = TRUE, spatial_network_name = 'kNN_network') Visualize top 2 ranked spatial genes per expression bin: set0 = ranktest[high_expr < 50][1:2]$feats set1 = ranktest[high_expr > 50 & high_expr < 100][1:2]$feats set2 = ranktest[high_expr > 100 & high_expr < 200][1:2]$feats set3 = ranktest[high_expr > 200 & high_expr < 400][1:2]$feats set4 = ranktest[high_expr > 400 & high_expr < 1000][1:2]$feats set5 = ranktest[high_expr > 1000][1:2]$feats spatFeatPlot2D(visiumHD, expression_values = 'scaled', feats = c(set0, set1, set2), gradient_style = "sequential", cell_color_gradient = c('blue', 'white', 'yellow', 'orange', 'red', 'darkred'), cow_n_col = 2, point_size = 1) Figure 9.16: Spat feature plot showing gene expression for the top 2 ranked spatial genes per expression bin (<50, >50 and >100) across the hex100 bin. spatFeatPlot2D(visiumHD, expression_values = 'scaled', feats = c(set3, set4, set5), gradient_style = "sequential", cell_color_gradient = c('blue', 'white', 'yellow', 'orange', 'red', 'darkred'), cow_n_col = 2, point_size = 1) Figure 9.17: Spat feature plot showing gene expression for the top 2 ranked spatial genes per expression bin (>200, >400 and >1000) across the hex100 bin. 9.5.2.2 Spatial co-expression modules ext_spatial_genes = ranktest[adj.p.value < 0.001]$feats spat_cor_netw_DT = detectSpatialCorFeats(visiumHD, method = 'network', spatial_network_name = 'kNN_network', subset_feats = ext_spatial_genes) # cluster spatial genes spat_cor_netw_DT = clusterSpatialCorFeats(spat_cor_netw_DT, name = 'spat_netw_clus', k = 16) # visualize clusters heatmSpatialCorFeats(visiumHD, spatCorObject = spat_cor_netw_DT, use_clus_name = 'spat_netw_clus', heatmap_legend_param = list(title = NULL)) Figure 9.18: Heatmap showing spatially correlated genes split into 16 clusters. # create metagene enrichment score for clusters cluster_genes_DT = showSpatialCorFeats(spat_cor_netw_DT, use_clus_name = 'spat_netw_clus', show_top_feats = 1) cluster_genes = cluster_genes_DT$clus; names(cluster_genes) = cluster_genes_DT$feat_ID visiumHD = createMetafeats(visiumHD, expression_values = 'normalized', feat_clusters = cluster_genes, name = 'cluster_metagene') showGiottoSpatEnrichments(visiumHD) spatCellPlot(visiumHD, spat_enr_names = 'cluster_metagene', gradient_style = "sequential", cell_color_gradient = c('blue', 'white', 'yellow', 'orange', 'red', 'darkred'), cell_annotation_values = as.character(c(1:4)), point_size = 1, cow_n_col = 2) Figure 9.19: Spat plot vizualizing metagenes (1-4) based on spatially correlated genes vizualized on the hex100 bin spatCellPlot(visiumHD, spat_enr_names = 'cluster_metagene', gradient_style = "sequential", cell_color_gradient = c('blue', 'white', 'yellow', 'orange', 'red', 'darkred'), cell_annotation_values = as.character(c(5:8)), point_size = 1, cow_n_col = 2) Figure 9.20: Spat plot vizualizing metagenes (5-8) based on spatially correlated genes vizualized on the hex100 bin spatCellPlot(visiumHD, spat_enr_names = 'cluster_metagene', gradient_style = "sequential", cell_color_gradient = c('blue', 'white', 'yellow', 'orange', 'red', 'darkred'), cell_annotation_values = as.character(c(9:12)), point_size = 1, cow_n_col = 2) Figure 9.21: Spat plot vizualizing metagenes (9-12) based on spatially correlated genes vizualized on the hex100 bin spatCellPlot(visiumHD, spat_enr_names = 'cluster_metagene', gradient_style = "sequential", cell_color_gradient = c('blue', 'white', 'yellow', 'orange', 'red', 'darkred'), cell_annotation_values = as.character(c(13:16)), point_size = 1, cow_n_col = 2) Figure 9.22: Spat plot vizualizing metagenes (13-16) based on spatially correlated genes vizualized on the hex100 bin 9.5.2.3 Plot spatial gene groups balanced_genes = getBalancedSpatCoexpressionFeats(spatCorObject = spat_cor_netw_DT, maximum = 5) selected_feats = names(balanced_genes) # give genes from same cluster same color distinct_colors = getDistinctColors(n = 20) names(distinct_colors) = 1:20 my_colors = distinct_colors[balanced_genes] names(my_colors) = names(balanced_genes) spatInSituPlotPoints(visiumHD, show_image = F, feats = list('rna' = selected_feats), feats_color_code = my_colors, show_legend = F, spat_unit = 'hex100', point_size = 0.20, show_polygon = FALSE, use_overlap = FALSE, polygon_feat_type = 'hex100', polygon_bg_color = NA, polygon_color = 'white', polygon_line_size = 0.01, expand_counts = TRUE, count_info_column = 'count', jitter = c(25,25)) Figure 9.23: Coloring individual features based on the spatially correlated gene clusters. 9.6 Hexbin 25 Goal is to create a higher resolution bin (hex25) and add to the Giotto object. We will aim to identify individual cell types and local neighborhood niches. 9.6.1 Subcellular workflow and projection functions filter and normalization workflow PCA projection visiumHD_subset = subsetGiottoLocs(gobject = visiumHD, x_min = 16000, x_max = 20000, y_min = 44250, y_max = 45500) spatInSituPlotPoints(visiumHD_subset, show_image = F, feats = NULL, show_legend = F, spat_unit = 'hex100', point_size = 0.5, show_polygon = TRUE, use_overlap = FALSE, polygon_feat_type = 'hex100', polygon_fill_as_factor = TRUE, polygon_fill = 'leiden_clus', polygon_color = 'black', polygon_line_size = 0.3) Figure 9.24: xxx spatInSituPlotPoints(visiumHD_subset, show_image = F, feats = list('rna' = selected_feats), feats_color_code = my_colors, show_legend = F, spat_unit = 'hex100', point_size = 0.40, show_polygon = TRUE, use_overlap = FALSE, polygon_feat_type = 'hex100', polygon_bg_color = NA, polygon_color = 'white', polygon_line_size = 0.05, jitter = c(25,25)) Figure 9.25: xxx hexbin25 <- tessellate(extent = ext(visiumHD_subset@feat_info$rna), shape = 'hexagon', shape_size = 25, name = 'hex25') visiumHD_subset = setPolygonInfo(gobject = visiumHD_subset, x = hexbin25, name = 'hex25', initialize = T) showGiottoSpatialInfo(visiumHD_subset) visiumHD_subset = addSpatialCentroidLocations(gobject = visiumHD_subset, poly_info = 'hex25') activeSpatUnit(visiumHD_subset) <- 'hex25' spatInSituPlotPoints(visiumHD_subset, show_image = F, feats = list('rna' = selected_feats), feats_color_code = my_colors, show_legend = F, spat_unit = 'hex25', point_size = 0.40, show_polygon = TRUE, use_overlap = FALSE, polygon_feat_type = 'hex25', polygon_bg_color = NA, polygon_color = 'white', polygon_line_size = 0.05, jitter = c(25,25)) Figure 9.26: xxx visiumHD_subset = calculateOverlap(visiumHD_subset, spatial_info = 'hex25', feat_info = 'rna') showGiottoSpatialInfo(visiumHD_subset) # convert overlap results to bin by gene matrix visiumHD_subset = overlapToMatrix(visiumHD_subset, poly_info = 'hex25', feat_info = 'rna', name = 'raw') visiumHD_subset <- filterGiotto(visiumHD_subset, expression_threshold = 1, feat_det_in_min_cells = 3, min_det_feats_per_cell = 5) activeSpatUnit(visiumHD_subset) # normalize visiumHD_subset <- normalizeGiotto(visiumHD_subset, scalefactor = 1000, verbose = T) # add statistics visiumHD_subset <- addStatistics(visiumHD_subset) feature_data = fDataDT(visiumHD_subset) visiumHD_subset <- calculateHVF(visiumHD_subset, zscore_threshold = 1) 9.6.2 projections 9.6.2.1 dimension reduction with projections n_25_percent <- round(length(spatIDs(visiumHD_subset, 'hex25')) * 0.25) # pca projection on subset visiumHD_subset <- runPCAprojection( gobject = visiumHD_subset, spat_unit = "hex25", feats_to_use = 'hvf', name = 'pca.projection', set_seed = TRUE, seed_number = 12345, random_subset = n_25_percent ) showGiottoDimRed(visiumHD_subset) plotPCA(visiumHD_subset, dim_reduction_name = 'pca.projection') Figure 9.27: xxx # umap projection on subset visiumHD_subset <- runUMAPprojection( gobject = visiumHD_subset, spat_unit = "hex25", dim_reduction_to_use = 'pca', dim_reduction_name = "pca.projection", dimensions_to_use = 1:10, name = "umap.projection", random_subset = n_25_percent, n_neighbors = 10, min_dist = 0.005, n_threads = 4 ) showGiottoDimRed(visiumHD_subset) # plot UMAP, coloring cells/points based on nr_feats plotUMAP(gobject = visiumHD_subset, point_size = 1, dim_reduction_name = 'umap.projection') Figure 9.28: xxx 9.6.2.2 clustering with projection # subset to smaller giotto object set.seed(1234) subset_IDs = sample(x = spatIDs(visiumHD_subset, 'hex25'), size = n_25_percent) temp_gobject = subsetGiotto( gobject = visiumHD_subset, spat_unit = 'hex25', cell_ids = subset_IDs ) # hierarchical clustering temp_gobject = doHclust(gobject = temp_gobject, spat_unit = 'hex25', k = 8, name = 'sub_hclust', dim_reduction_to_use = 'pca', dim_reduction_name = 'pca.projection', dimensions_to_use = 1:10) # show umap dimPlot2D( gobject = temp_gobject, point_size = 2.5, spat_unit = 'hex25', dim_reduction_to_use = 'umap', dim_reduction_name = 'umap.projection', cell_color = 'sub_hclust' ) Figure 9.29: xxx # project clusterings back to full dataset visiumHD_subset <- doClusterProjection( target_gobject = visiumHD_subset, source_gobject = temp_gobject, spat_unit = "hex25", source_cluster_labels = "sub_hclust", reduction_method = 'pca', reduction_name = 'pca.projection', prob = FALSE, knn_k = 5, dimensions_to_use = 1:10 ) pDataDT(visiumHD_subset) dimPlot2D( gobject = visiumHD_subset, point_size = 1.5, spat_unit = 'hex25', dim_reduction_to_use = 'umap', dim_reduction_name = 'umap.projection', cell_color = 'knn_labels' ) Figure 9.30: xxx spatInSituPlotPoints(visiumHD_subset, show_image = F, feats = NULL, show_legend = F, spat_unit = 'hex25', point_size = 0.5, show_polygon = TRUE, use_overlap = FALSE, polygon_feat_type = 'hex25', polygon_fill_as_factor = TRUE, polygon_fill = 'knn_labels', polygon_color = 'black', polygon_line_size = 0.3) Figure 9.31: xxx 9.6.3 niche clusters visiumHD_subset = createSpatialNetwork(visiumHD_subset, name = 'kNN_network', spat_unit = 'hex25', method = 'kNN', k = 6) pDataDT(visiumHD_subset) visiumHD_subset = calculateSpatCellMetadataProportions(gobject = visiumHD_subset, spat_unit = 'hex25', feat_type = 'rna', metadata_column = 'knn_labels', spat_network = 'kNN_network') prop_table = getSpatialEnrichment(visiumHD_subset, name = 'proportion', output = 'data.table') prop_matrix = GiottoUtils:::dt_to_matrix(prop_table) set.seed(1234) prop_kmeans = kmeans(x = prop_matrix, centers = 10, iter.max = 1000, nstart = 100) prop_kmeansDT = data.table::data.table(cell_ID = names(prop_kmeans$cluster), niche = prop_kmeans$cluster) visiumHD_subset = addCellMetadata(visiumHD_subset, new_metadata = prop_kmeansDT, by_column = T, column_cell_ID = 'cell_ID') pDataDT(visiumHD_subset) spatInSituPlotPoints(visiumHD_subset, show_image = F, feats = NULL, show_legend = F, spat_unit = 'hex25', point_size = 0.5, show_polygon = TRUE, use_overlap = FALSE, polygon_feat_type = 'hex25', polygon_fill_as_factor = TRUE, polygon_fill = 'niche', polygon_color = 'black', polygon_line_size = 0.3) Figure 9.32: xxx 9.7 Database backend - Work in progress, but coming soon! Memory problems: - data ingestion - spatial operations - matrix operations - matrix and spatial geometry object sizes "],["404.html", "Page not found", " Page not found The page you requested cannot be found (perhaps it was moved or renamed). You may want to try searching to find the page's new location, or use the table of contents to find the page you are looking for. "]] diff --git a/docs/visium-hd.html b/docs/visium-hd.html index fb79709..bec6ab7 100644 --- a/docs/visium-hd.html +++ b/docs/visium-hd.html @@ -325,6 +325,8 @@
  • 9.6 Hexbin 25
  • 9.7 Database backend - Work in progress, but coming soon!
    • @@ -632,8 +634,8 @@

    9.3.2 Mini Visium HD dataset# download the mini dataset and untar options("timeout" = Inf) download.file( - url = "https://zenodo.org/records/13212855/files/workshop_VisiumHD.zip?download=1", - destfile = file.path(save_dir, "workshop_visiumHD.zip") + url = "https://zenodo.org/records/13226158/files/workshop_VisiumHD.zip?download=1", + destfile = file.path(save_dir, "workshop_visiumHD.zip") ) untar(tarfile = file.path(save_dir, "workshop_visiumHD.zip"), exdir = data_path) @@ -648,13 +650,13 @@

    9.3.4 Read in data manuallyHowever, for this tutorial we will illustrate how to create your own Giotto object in a step-by-step manner, which can also be applied to other similar technologies as discussed in the Objective section.

    9.3.4.1 Raw expression data

    -
    expression_path <- file.path(data_path, 'workshop_visiumHD/Human_Colorectal_Cancer_workshop/square_002um/raw_feature_bc_matrix')
+
    expression_path <- file.path(data_path, '/Human_Colorectal_Cancer_workshop/square_002um/raw_feature_bc_matrix')
 expr_results <- get10Xmatrix(path_to_data = expression_path, 
                              gene_column_index = 1)

    9.3.4.2 Tissue positions data

    -
    tissue_positions_path <- file.path(data_path, 'workshop_visiumHD/Human_Colorectal_Cancer_workshop/square_002um/spatial/tissue_positions.parquet')
+
    tissue_positions_path <- file.path(data_path, '/Human_Colorectal_Cancer_workshop/square_002um/spatial/tissue_positions.parquet')
 tissue_positions <- data.table::as.data.table(arrow::read_parquet(tissue_positions_path))
    @@ -728,9 +730,7 @@

    9.4.3 combine Giotto points and p visiumHD = addSpatialCentroidLocations(gobject = visiumHD, poly_info = 'hex400')

    -

    Visualize the Giotto object. Make sure to set expand_counts = TRUE to expand the counts -column. Each spatial bin can have multiple transcripts/UMIs. This is different -compared to in situ technologies like seqFISH, MERFISH, Nanostring CosMx or Xenium.

    +

    Visualize the Giotto object. Make sure to set expand_counts = TRUE to expand the counts column. Each spatial bin can have multiple transcripts/UMIs. This is different compared to in situ technologies like seqFISH, MERFISH, Nanostring CosMx or Xenium.

    feature_data = fDataDT(visiumHD)
 
 spatInSituPlotPoints(visiumHD,
@@ -754,8 +754,7 @@ 
    9.4.3 combine Giotto points and p
 Figure 9.7: Overlap of gene expression with the hex400 polygons. Each dot represents a single gene. Jitter used to better vizualize individual genes.
    -

    You can set plot_method = scattermore or scattermost to convert high-resolution -images to low(er) resolution rasterized images. It’s usually faster and will save on disk space.

    +

    You can set plot_method = scattermore or scattermost to convert high-resolution images to low(er) resolution rasterized images. It’s usually faster and will save on disk space.

    spatInSituPlotPoints(visiumHD,
                      show_image = F,
                      feats = list('rna' = feature_data$feat_ID[10:20]),
@@ -1210,8 +1209,7 @@ 
    9.5.2.3 Plot spatial gene groups<

    9.6 Hexbin 25

    -

    Goal is to create a higher resolution bin (hex25) and add to the Giotto object. -We will aim to identify individual cell types and local neighborhood niches.

    +

    Goal is to create a higher resolution bin (hex25) and add to the Giotto object. We will aim to identify individual cell types and local neighborhood niches.

    9.6.1 Subcellular workflow and projection functions

      @@ -1223,44 +1221,308 @@

      9.6.1 Subcellular workflow and pr x_max = 20000, y_min = 44250, y_max = 45500)

    -
    hexbin25 <- tessellate(extent = ext(visiumHD_subset@feat_info$rna), 
-                       shape = 'hexagon', 
-                       shape_size = 25, 
-                       name = 'hex25') 
-
-visiumHD_subset = setPolygonInfo(gobject = visiumHD_subset,
-                                 x = hexbin25,
-                                 name = 'hex25',
-                                 initialize = T)
-
-visiumHD_subset = addSpatialCentroidLocations(gobject = visiumHD_subset,
-                                              poly_info = 'hex25')
-
-activeSpatUnit(visiumHD_subset) <- 'hex25'
    -
    visiumHD_subset = calculateOverlap(visiumHD_subset,
-                                   spatial_info = 'hex25',
-                                   feat_info = 'rna')
-
-# convert overlap results to bin by gene matrix
-visiumHD_subset = overlapToMatrix(visiumHD_subset,
-                                  poly_info = 'hex25',
-                                  feat_info = 'rna',
-                                  name = 'raw')
-
-visiumHD_subset <- filterGiotto(visiumHD_subset,
-                                expression_threshold = 1,
-                                feat_det_in_min_cells = 25,
-                                min_det_feats_per_cell = 25)
-
-# normalize
-visiumHD_subset <- normalizeGiotto(visiumHD_subset, scalefactor = 1000, verbose = T)
-
-# add statistics
-visiumHD_subset <- addStatistics(visiumHD_subset)
    +
    spatInSituPlotPoints(visiumHD_subset,
+                     show_image = F,
+                     feats = NULL,
+                     show_legend = F,
+                     spat_unit = 'hex100',
+                     point_size = 0.5,
+                     show_polygon = TRUE,
+                     use_overlap = FALSE,
+                     polygon_feat_type = 'hex100',
+                     polygon_fill_as_factor = TRUE,
+                     polygon_fill = 'leiden_clus',
+                     polygon_color = 'black',
+                     polygon_line_size = 0.3)
    +
    +xxx +

    +Figure 9.24: xxx +

    +
    +
    spatInSituPlotPoints(visiumHD_subset,
+                     show_image = F,
+                     feats = list('rna' = selected_feats), 
+                     feats_color_code = my_colors,
+                     show_legend = F,
+                     spat_unit = 'hex100',
+                     point_size = 0.40,
+                     show_polygon = TRUE,
+                     use_overlap = FALSE,
+                     polygon_feat_type = 'hex100',
+                     polygon_bg_color = NA,
+                     polygon_color = 'white',
+                     polygon_line_size = 0.05,
+                     jitter = c(25,25))
    +
    +xxx +

    +Figure 9.25: xxx +

    +
    +
    hexbin25 <- tessellate(extent = ext(visiumHD_subset@feat_info$rna), 
+                       shape = 'hexagon', 
+                       shape_size = 25, 
+                       name = 'hex25') 
+
+visiumHD_subset = setPolygonInfo(gobject = visiumHD_subset,
+                                 x = hexbin25,
+                                 name = 'hex25',
+                                 initialize = T)
+
+showGiottoSpatialInfo(visiumHD_subset)
+
+visiumHD_subset = addSpatialCentroidLocations(gobject = visiumHD_subset,
+                                              poly_info = 'hex25')
+
+activeSpatUnit(visiumHD_subset) <- 'hex25'
+
+spatInSituPlotPoints(visiumHD_subset,
+                     show_image = F,
+                     feats = list('rna' = selected_feats), 
+                     feats_color_code = my_colors,
+                     show_legend = F,
+                     spat_unit = 'hex25',
+                     point_size = 0.40,
+                     show_polygon = TRUE,
+                     use_overlap = FALSE,
+                     polygon_feat_type = 'hex25',
+                     polygon_bg_color = NA,
+                     polygon_color = 'white',
+                     polygon_line_size = 0.05,
+                     jitter = c(25,25))
    +
    +xxx +

    +Figure 9.26: xxx +

    +
    +
    visiumHD_subset = calculateOverlap(visiumHD_subset,
+                                   spatial_info = 'hex25',
+                                   feat_info = 'rna')
+
+showGiottoSpatialInfo(visiumHD_subset)
+
+
+# convert overlap results to bin by gene matrix
+visiumHD_subset = overlapToMatrix(visiumHD_subset,
+                                  poly_info = 'hex25',
+                                  feat_info = 'rna',
+                                  name = 'raw')
+
+visiumHD_subset <- filterGiotto(visiumHD_subset,
+                                expression_threshold = 1,
+                                feat_det_in_min_cells = 3,
+                                min_det_feats_per_cell = 5)
+
+activeSpatUnit(visiumHD_subset)
+
+
+# normalize
+visiumHD_subset <- normalizeGiotto(visiumHD_subset, scalefactor = 1000, verbose = T)
+
+# add statistics
+visiumHD_subset <- addStatistics(visiumHD_subset)
+
+feature_data = fDataDT(visiumHD_subset)
+
+visiumHD_subset <- calculateHVF(visiumHD_subset, zscore_threshold = 1)
    +
    +
    +

    9.6.2 projections

    +
    +

    9.6.2.1 dimension reduction with projections

    +
    n_25_percent <- round(length(spatIDs(visiumHD_subset, 'hex25')) * 0.25)
+
+# pca projection on subset
+visiumHD_subset <- runPCAprojection(
+  gobject = visiumHD_subset,
+  spat_unit = "hex25",
+  feats_to_use = 'hvf',
+  name = 'pca.projection',
+  set_seed = TRUE,
+  seed_number = 12345,
+  random_subset = n_25_percent
+)
+
+showGiottoDimRed(visiumHD_subset)
+plotPCA(visiumHD_subset, dim_reduction_name = 'pca.projection')
    +
    +xxx +

    +Figure 9.27: xxx +

    +
    +
    # umap projection on subset
+visiumHD_subset <- runUMAPprojection(
+  gobject = visiumHD_subset,
+  spat_unit = "hex25",
+  dim_reduction_to_use = 'pca',
+  dim_reduction_name = "pca.projection",
+  dimensions_to_use = 1:10,
+  name = "umap.projection",
+  random_subset = n_25_percent, 
+  n_neighbors = 10,
+  min_dist = 0.005,
+  n_threads = 4
+)
+
+showGiottoDimRed(visiumHD_subset)
+
+# plot UMAP, coloring cells/points based on nr_feats
+plotUMAP(gobject = visiumHD_subset,
+         point_size = 1, 
+         dim_reduction_name = 'umap.projection')
    +
    +xxx +

    +Figure 9.28: xxx +

    +
    +
    +
    +

    9.6.2.2 clustering with projection

    +
    # subset to smaller giotto object
+set.seed(1234)
+subset_IDs = sample(x = spatIDs(visiumHD_subset, 'hex25'), size = n_25_percent)
+
+temp_gobject = subsetGiotto(
+  gobject = visiumHD_subset,
+  spat_unit = 'hex25',
+  cell_ids = subset_IDs
+)
+
+
+# hierarchical clustering
+temp_gobject = doHclust(gobject = temp_gobject, 
+                        spat_unit = 'hex25', 
+                        k = 8, name = 'sub_hclust', 
+                        dim_reduction_to_use = 'pca', 
+                        dim_reduction_name = 'pca.projection', 
+                        dimensions_to_use = 1:10)
+
+
+# show umap
+dimPlot2D(
+  gobject = temp_gobject,
+  point_size = 2.5,
+  spat_unit = 'hex25',
+  dim_reduction_to_use = 'umap',
+  dim_reduction_name = 'umap.projection',
+  cell_color = 'sub_hclust'
+)
    +
    +xxx +

    +Figure 9.29: xxx +

    +
    +
    # project clusterings back to full dataset
+visiumHD_subset <- doClusterProjection(
+  target_gobject = visiumHD_subset,
+  source_gobject = temp_gobject,
+  spat_unit = "hex25",
+  source_cluster_labels = "sub_hclust",
+  reduction_method = 'pca',
+  reduction_name = 'pca.projection',
+  prob = FALSE,
+  knn_k = 5,
+  dimensions_to_use = 1:10
+)
+
+pDataDT(visiumHD_subset)
+
+dimPlot2D(
+  gobject = visiumHD_subset,
+  point_size = 1.5,
+  spat_unit = 'hex25',
+  dim_reduction_to_use = 'umap',
+  dim_reduction_name = 'umap.projection',
+  cell_color = 'knn_labels'
+)
    +
    +xxx +

    +Figure 9.30: xxx +

    +
    +
    spatInSituPlotPoints(visiumHD_subset,
+                     show_image = F,
+                     feats = NULL,
+                     show_legend = F,
+                     spat_unit = 'hex25',
+                     point_size = 0.5,
+                     show_polygon = TRUE,
+                     use_overlap = FALSE,
+                     polygon_feat_type = 'hex25',
+                     polygon_fill_as_factor = TRUE,
+                     polygon_fill = 'knn_labels',
+                     polygon_color = 'black',
+                     polygon_line_size = 0.3)
    +
    +xxx +

    +Figure 9.31: xxx +

    +
    +
    +
    +
    +

    9.6.3 niche clusters

    +
    visiumHD_subset = createSpatialNetwork(visiumHD_subset,
+                                name = 'kNN_network',
+                                spat_unit = 'hex25', 
+                                method = 'kNN',
+                                k = 6)
+
+pDataDT(visiumHD_subset)
+visiumHD_subset = calculateSpatCellMetadataProportions(gobject = visiumHD_subset,
+                                                       spat_unit = 'hex25',
+                                                       feat_type = 'rna',
+                                                       metadata_column = 'knn_labels',
+                                                       spat_network = 'kNN_network')
+
+prop_table = getSpatialEnrichment(visiumHD_subset, name = 'proportion', output = 'data.table')
+prop_matrix = GiottoUtils:::dt_to_matrix(prop_table)
+
+set.seed(1234)
+prop_kmeans = kmeans(x = prop_matrix, centers = 10, iter.max = 1000, nstart = 100)
+prop_kmeansDT = data.table::data.table(cell_ID = names(prop_kmeans$cluster), niche = prop_kmeans$cluster)
+
+visiumHD_subset = addCellMetadata(visiumHD_subset, 
+                                  new_metadata = prop_kmeansDT, 
+                                  by_column = T, 
+                                  column_cell_ID = 'cell_ID')
+pDataDT(visiumHD_subset)
+
+spatInSituPlotPoints(visiumHD_subset,
+                     show_image = F,
+                     feats = NULL,
+                     show_legend = F,
+                     spat_unit = 'hex25',
+                     point_size = 0.5,
+                     show_polygon = TRUE,
+                     use_overlap = FALSE,
+                     polygon_feat_type = 'hex25',
+                     polygon_fill_as_factor = TRUE,
+                     polygon_fill = 'niche',
+                     polygon_color = 'black',
+                     polygon_line_size = 0.3)
    +
    +xxx +

    +Figure 9.32: xxx +

    +

    9.7 Database backend - Work in progress, but coming soon!

    +

    Memory problems:
    +- data ingestion
    +- spatial operations
    +- matrix operations
    +- matrix and spatial geometry object sizes

    diff --git a/img/02_session2/0-spatInSituPlotPoints.png b/img/02_session2/0-spatInSituPlotPoints.png index 83a6f19..4832d47 100644 Binary files a/img/02_session2/0-spatInSituPlotPoints.png and b/img/02_session2/0-spatInSituPlotPoints.png differ diff --git a/img/02_session2/1-spatInSituPlotPoints.png b/img/02_session2/1-spatInSituPlotPoints.png index 58ddb68..a769766 100644 Binary files a/img/02_session2/1-spatInSituPlotPoints.png and b/img/02_session2/1-spatInSituPlotPoints.png differ diff --git a/img/02_session2/12-UMAP.png b/img/02_session2/12-UMAP.png index 3e6e878..bbbb2ae 100644 Binary files a/img/02_session2/12-UMAP.png and b/img/02_session2/12-UMAP.png differ diff --git a/img/02_session2/21-spatInSituPlotPoints.png b/img/02_session2/21-spatInSituPlotPoints.png index 47dcfb8..96f4ac4 100644 Binary files a/img/02_session2/21-spatInSituPlotPoints.png and b/img/02_session2/21-spatInSituPlotPoints.png differ diff --git a/img/02_session2/22-spatInSituPlotPoints.png b/img/02_session2/22-spatInSituPlotPoints.png new file mode 100644 index 0000000..490010a Binary files /dev/null and b/img/02_session2/22-spatInSituPlotPoints.png differ diff --git a/img/02_session2/23-spatInSituPlotPoints.png b/img/02_session2/23-spatInSituPlotPoints.png new file mode 100644 index 0000000..c77c516 Binary files /dev/null and b/img/02_session2/23-spatInSituPlotPoints.png differ diff --git a/img/02_session2/24-spatInSituPlotPoints.png b/img/02_session2/24-spatInSituPlotPoints.png new file mode 100644 index 0000000..fb4eb71 Binary files /dev/null and b/img/02_session2/24-spatInSituPlotPoints.png differ diff --git a/img/02_session2/25-PCA.png b/img/02_session2/25-PCA.png new file mode 100644 index 0000000..0620e77 Binary files /dev/null and b/img/02_session2/25-PCA.png differ diff --git a/img/02_session2/26-UMAP.png b/img/02_session2/26-UMAP.png new file mode 100644 index 0000000..7b5f766 Binary files /dev/null and b/img/02_session2/26-UMAP.png differ diff --git a/img/02_session2/27-dimPlot2D.png b/img/02_session2/27-dimPlot2D.png new file mode 100644 index 0000000..0c65032 Binary files /dev/null and b/img/02_session2/27-dimPlot2D.png differ diff --git a/img/02_session2/28-dimPlot2D.png b/img/02_session2/28-dimPlot2D.png new file mode 100644 index 0000000..0764066 Binary files /dev/null and b/img/02_session2/28-dimPlot2D.png differ diff --git a/img/02_session2/29-spatInSituPlotPoints.png b/img/02_session2/29-spatInSituPlotPoints.png new file mode 100644 index 0000000..c65d5a2 Binary files /dev/null and b/img/02_session2/29-spatInSituPlotPoints.png differ diff --git a/img/02_session2/30-spatInSituPlotPoints.png b/img/02_session2/30-spatInSituPlotPoints.png new file mode 100644 index 0000000..9c9f871 Binary files /dev/null and b/img/02_session2/30-spatInSituPlotPoints.png differ diff --git a/img/02_session2/7-UMAP.png b/img/02_session2/7-UMAP.png index 2c2cd83..95215d1 100644 Binary files a/img/02_session2/7-UMAP.png and b/img/02_session2/7-UMAP.png differ diff --git a/img/02_session2/9-spatInSituPlotPoints.png b/img/02_session2/9-spatInSituPlotPoints.png index cb11008..723b9d1 100644 Binary files a/img/02_session2/9-spatInSituPlotPoints.png and b/img/02_session2/9-spatInSituPlotPoints.png differ