How to properly subset by FOV

[rbutleriii]

I have a CosMx slide that looks at multiple tissue sections per slide (different experimental groups). I want to split the slide by each brain slice, which itself has multiple fovs:

This can be easily accomplished on the cell_ID side with expression metadata:

# subset by group
cells = pDataDT(gobj)[group == grp, cell_ID]
gobj = subsetGiotto(gobj, cell_ids=cells)

However the images, LargeImages and the global subcellular features remain. So even though I can plot the single slice without images

(or even with images if just looking at cell polygons)

once I plot the features, all the other fovs are displayed

Is there a way to drop images and/or features points by fov?

[rbutleriii]

So, trying to make some progress on a workaround, and it looks like it is easy enough to drop the extra images:

# dropping images
gobj@images = gobj@images[1]
gobj@largeImages = gobj@largeImages[1:3]

However, the issue is the features/polygons. Looks like I am then trying to filter spatial_info, spatial_locs, polygon_info, etc all one at a time. Does subsetGiottoLocs handle all the requisite feature/polygon/image processing as well? If so, I would just need to grab the xmin/xmax/ymin/ymax for each image boundary...

A thought just occurred, if I subset the expression matrix, and the features/cell_polygons remain, are the downstream analyses (spatial network, cell-cell interaction, etc) going to include any features outside of my subset expression matrix?

[RubD]

@rbutleriii I believe that @jiajic has been working on a tutorial on the different options of subsetting spatial data and we hope to put that online in the near future. For the images we currently don't really have a solution and your workaround setting those to NULL is probably a good idea for now.

[rbutleriii]

Got it, I have been reworking my framing of the prep to avoid needing to subset the fovs, which is easy enough by making separate giotto objects for each brain slice.

Ultimately, when subsetting comes online I would end up wanting to subset by region and/or subregion (in this case cortex vs hippocampus, and then CA1 vs CA23 vs DG), which is better accomplished using an interactive selection session to outline those regions. If the subcellular feature data also obeyed those polygon bounds it would make for a good image (as in the tutorials, with points and cell polygons). I also took a look at the latest GiottoClass version (0.1.3) of subsetGiotto, and I think there looks like there are options to make that work, I am just not sure how to set them, something like:

library(Giotto)
library(GiottoData)
library(data.table)

x = loadGiottoMini("cosmx")

# Filtering & normalization
x %>%
  filterGiotto(
    feat_type = 'rna',
    feat_det_in_min_cells = 5,
    min_det_feats_per_cell = 5
  ) %>%
  normalizeGiotto(
    feat_type = 'rna'
  ) %>%
  normalizeGiotto(
    feat_type = 'neg_probe',
    library_size_norm = FALSE
  ) %>%
  normalizeGiotto(
    feat_type = 'rna',
    scalefactor = 5000,
    norm_methods = 'pearson_resid',
    update_slot = 'pearson'
  ) -> x

# Plot spatial cells/spots with tissue on the background
x %>%
  spatPlot2D(
    point_size = 1.5,
    return_plot = TRUE,
    show_plot = FALSE
  ) %>%
  plotInteractivePolygons() -> my_polygon_coordinates

# > dput(my_polygon_coordinates)
my_polygon_coordinates <- structure(
  list(x = c(5367.32407979834, 5369.22771954729, 5989.81427770447, 
  6025.98343293449, 5822.29397979701, 6168.75641410562, 6115.45450113507, 
  6625.62995385324, 6646.56999109168), y = c(-156663.925097282, 
  -157155.064152511, -157231.209742469, -156633.466861299, -156603.008625316, 
  -156757.203444981, -157195.040587239, -157235.017021967, -156877.132749165
  ), name = c("polygon 1", "polygon 1", "polygon 1", "polygon 1", 
  "polygon 1", "polygon 2", "polygon 2", "polygon 2", "polygon 2"
  )), row.names = c(NA, -9L), class = c("data.table", "data.frame")
)

# We must transform the data.table or data.frame with coordinates into a Giotto polygon object
my_giotto_polygons <- createGiottoPolygonsFromDfr(
  my_polygon_coordinates, 
  name = 'selections'
)

# Then, add the polygons to the Giotto object
x %>% 
  addGiottoPolygons(gpolygons = list(my_giotto_polygons)) %>%
  addPolygonCells(polygon_name = 'selections') -> x

# plot some features everywhere
hvfs = fDataDT(x)[order(-var)][hvf == 'yes', feat_ID]
spatInSituPlotPoints(x,
  feats = list(hvfs[1:min(5, length(hvfs))]),
  point_size = 1,
  polygon_line_size = 0.1,
  show_polygon = TRUE,
  polygon_color = 'white',
  show_image = FALSE,
  # largeImage_name = names(x@largeImages),
  show_legend = TRUE
)

# plot some features only in the bounds of the polygon 1 region 
y = subsetGiotto(x, spat_unit_fsub="polygon 1")
spatInSituPlotPoints(y,
  feats = list(hvfs[1:min(5, length(hvfs))]),
  point_size = 1,
  polygon_line_size = 0.1,
  show_polygon = TRUE,
  polygon_color = 'white',
  show_image = FALSE,
  # largeImage_name = names(x@largeImages),
  show_legend = TRUE
)