library(Giotto)

# Ensure Giotto can access a python env
genv_exists <- suppressMessages(checkGiottoEnvironment())
print(genv_exists)

python_path <- file.path(Sys.getenv("HOME"), "miniconda3", "envs", "giotto_env", "bin", "python")
Sys.setenv(RETICULATE_PYTHON = python_path)

invisible(capture.output(suppressMessages(set_giotto_python_path(python_path = python_path))))

# library(scplotter)
devtools::load_all()

Loading required package: GiottoClass

Newer devel version of GiottoClass on GitHub: 0.4.8

Giotto Suite 4.2.1

[1] TRUE

ℹ Loading scplotter

library(Giotto)

## Set instructions
results_folder <- "data/Nanostring_CosMx_Subcellular_Lung_Cancer.results"

instrs <- createGiottoInstructions(
    save_dir = results_folder,
    save_plot = FALSE,
    show_plot = TRUE,
    return_plot = TRUE,
    python_path = python_path
)

data_path <- "data/Nanostring_CosMx_Subcellular_Lung_Cancer/Lung12/Lung12-Flat_files_and_images/"

## create giotto cosmx object
cosmx <- suppressWarnings(createGiottoCosMxObject(
    cosmx_dir = data_path,
    version = "legacy", # set to this for legacy NSCLC dataset
    FOVs = seq_len(28), # defaults to all FOVs
    load_expression = FALSE, # defaults to FALSE (see feature aggregation step below)
    instructions = instrs
))

force(cosmx)

python already initialized in this session
 active environment : 'giotto_env'
 python version : 3.10

python already initialized in this session
 active environment : 'giotto_env'
 python version : 3.10

loading feature detections...

Setting feature info [rna]

Setting feature info [negprobes]

loading composite images...

loading overlay images...

loading cell metadata...

An object of class giotto 
>Active spat_unit:  cell 
>Active feat_type:  rna 
dimensions    : 960, 73997 (features, cells)
[SUBCELLULAR INFO]
polygons      : cell 
features      : rna negprobes 
[AGGREGATE INFO]
spatial locations ----------------
  [cell] raw
attached images ------------------
images      : 56 items...


Use objHistory() to see steps and params used

options(repr.plot.width = 10, repr.plot.height = 10)

SpatDimPlot(
    cosmx,
    label = TRUE,
    group_by = "fov",
    points_size = 0.1,
    image = "black"
)

options(repr.plot.width = 11, repr.plot.height = 10)
# devtools::load_all()

SpatDimPlot(
    cosmx,
    image = TRUE,  # use all fovs
    features = c(
        "MMP2", "VEGFA", "IGF1R",
        "MKI67", "EPCAM", "KRT8"
    ),
    spat_unit = "cell",
    points_size = 0.01,
    points_palette = "Vivid"
)

options(repr.plot.width = 11, repr.plot.height = 10)
# Visualize Cell Centroids

SpatDimPlot(
    cosmx,
    image = TRUE,  # use all fovs
    points_size = 0.1,
    points_alpha = 0.8,
    points_color_by = "white",
    title = "centroids"
)

# Aggregate subcellular features
# Find the feature points overlapped by polygons. This overlap information is then
# returned to the relevant giottoPolygon object"s overlaps slot.
cosmx <- calculateOverlapRaster(cosmx, feat_info = "rna")
cosmx <- calculateOverlapRaster(cosmx, feat_info = "negprobes")

# Convert the overlap information into a cell by feature expression matrix which
# is then stored in the Giotto object"s expression slot
cosmx <- overlapToMatrix(cosmx, feat_info = "rna")
cosmx <- overlapToMatrix(cosmx, feat_info = "negprobes")

cosmx <- filterGiotto(cosmx,
    feat_type = "rna",
    expression_threshold = 1,
    feat_det_in_min_cells = 5,
    min_det_feats_per_cell = 5
)

# normalize
# standard method of normalization (log normalization based)
cosmx <- normalizeGiotto(cosmx,
    feat_type = "rna",
    norm_methods = "standard"
)

# new normalization method based on pearson correlations (Lause/Kobak et al. 2021)
# this normalized matrix is given the name "pearson" using the update_slot param
cosmx <- normalizeGiotto(cosmx,
    feat_type = "rna",
    scalefactor = 5000,
    norm_methods = "pearson_resid",
    name = "pearson"
)
# add statistics based on log normalized values
cosmx <- addStatistics(cosmx,
    expression_values = "normalized",
    feat_type = "rna"
)

# add statistics for the raw negative probe values
cosmx <- addStatistics(cosmx,
    expression_values = "raw",
    feat_type = "negprobes"
)

# View cellular data (default is feat = "rna")
showGiottoCellMetadata(cosmx)

# View feature data
showGiottoFeatMetadata(cosmx)

1. convert polygon to raster

2. overlap raster and points

3. add polygon information

4. add points information

5. create overlap polygon
 information

1. convert polygon to raster

2. overlap raster and points

3. add polygon information

4. add points information

5. create overlap polygon
 information

completed 1: preparation

completed 2: subset expression data

completed 3: subset spatial locations

completed 4: subset cell metadata

completed 5: subset feature metadata

completed 6: subset spatial network(s)

completed 7: subsetted dimension reductions

completed 8: subsetted nearest network(s)

completed 9: subsetted spatial enrichment results

for cell

--> cell found back in polygon layer: cell

completed 10: subsetted spatial information data

subset feature info: rna

completed 11: subsetted spatial feature data

Feature type:  rna 
Number of cells removed:  775  out of  73997 
Number of feats removed:  0  out of  960

first scale feats and then cells

Setting expression [cell][rna] normalized

Setting expression [cell][rna] scaled

using 'Lause/Kobak' method to normalize count matrix If used in
        published research, please cite:
        Jan Lause, Philipp Berens, Dmitry Kobak (2020).
        'Analytic Pearson residuals for normalization of single-cell RNA-seq UMI
        data' 

Setting expression [cell][rna] pearson

calculating statistics for "normalized" expression

calculating statistics for "raw" expression


These column names were already used: area
 and will be overwritten

└──Spatial unit "cell"
   ├──Feature type "rna"
   │     An object of class cellMetaObj 
   │     spat_unit : "cell"
   │     feat_type : "rna"
   │     provenance: cell 
   │     dimensions: 4 20 
   │     
   │        cell_ID  Area AspectRatio Width Height Mean.MembraneStain Max.MembraneStain
   │         <char> <int>       <num> <int>  <int>              <int>             <int>
   │     1: c_1_1_1  3411        1.98    93     47               7127             16033
   │     2: c_1_1_2  2993        1.77    78     44              13547             25859
   │     3: c_1_1_3  5590        1.75   112     64              10353             22954
   │        Mean.PanCK Max.PanCK Mean.CD45 Max.CD45 Mean.CD3 Max.CD3 Mean.DAPI Max.DAPI
   │             <int>     <int>     <int>    <int>    <int>   <int>     <int>    <int>
   │     1:       3000     17172        78      588       80     610     12242    47623
   │     2:      12494     23204       102      410        0      85         0        0
   │     3:       7951     21297        79      495        0     133       611    23078
   │           fov nr_feats perc_feats total_expr  area
   │        <char>    <int>      <num>      <num> <num>
   │     1:      1       82   8.541667   493.4018  3411
   │     2:      1       52   5.416667   347.6662  2993
   │     3:      1       75   7.812500   467.2596  5590
   │  
   └──Feature type "negprobes"
         An object of class cellMetaObj 
         spat_unit : "cell"
         feat_type : "negprobes"
         provenance: cell 
         dimensions: 4 4 
         
            cell_ID nr_feats perc_feats total_expr
             <char>    <int>      <num>      <num>
         1: c_1_1_1        1          5          1
         2: c_1_1_2        0          0          0
         3: c_1_1_3        0          0          0
      
└──Spatial unit "cell"
   ├──Feature type "rna"
   │     An object of class featMetaObj 
   │     spat_unit : "cell"
   │     feat_type : "rna"
   │     provenance: cell 
   │     
   │        feat_ID nr_cells perc_cells total_expr mean_expr mean_expr_det
   │         <char>    <int>      <num>      <num>     <num>         <num>
   │     1:    IL7R     6700   9.150255   31279.27 0.4271840      4.668548
   │     2:  SEC61G     9009  12.303679   40746.90 0.5564844      4.522910
   │     3:   IGHA1    19286  26.339078  107382.91 1.4665389      5.567920
   │  
   └──Feature type "negprobes"
         An object of class featMetaObj 
         spat_unit : "cell"
         feat_type : "negprobes"
         provenance: cell 
         
             feat_ID nr_cells perc_cells total_expr  mean_expr mean_expr_det
              <char>    <int>      <num>      <num>      <num>         <num>
         1: NegPrb15     3257   4.448117       3674 0.05017618      1.128032
         2: NegPrb18     4178   5.705935       4799 0.06554041      1.148636
         3:  NegPrb7     3866   5.279834       4520 0.06173008      1.169167

options(repr.plot.width = 11, repr.plot.height = 10)

SpatFeaturePlot(
    cosmx,
    image = "black",
    feature = "total_expr",
    points_size = 0.4,
    points_palette = "BluYl",
    points_palette_reverse = TRUE
)

options(repr.plot.width = 11, repr.plot.height = 10)

SpatFeaturePlot(
    cosmx,
    image = "black",
    points = FALSE,
    shapes_fill_by = "total_expr",  # indicating shapes = TRUE
    shapes_border_size = 0.05
)

cosmx <- calculateHVF(cosmx,
    method = "var_p_resid",
    var_number = 500, # leaving as NULL results in too few genes to use
    save_plot = TRUE,
    save_param = list(
        save_name = "10_hvf"
    )
)
# If you get an Error related to future.apply, please modify the maximum size
# of global variables by running: options(future.globals.maxSize = 1e10)

# print HVFs
gene_metadata <- fDataDT(cosmx)
gene_metadata[hvf == "yes", feat_ID]

cosmx <- runPCA(cosmx,
    scale_unit = FALSE,
    center = FALSE,
    expression_values = "pearson"
)

cosmx <- runtSNE(cosmx,
    dimensions_to_use = 1:10
)

cosmx <- runUMAP(cosmx,
    dimensions_to_use = 1:10
)

cosmx <- createNearestNetwork(cosmx,
    dimensions_to_use = 1:10,
    k = 10
)

cosmx <- doLeidenCluster(cosmx,
    resolution = 0.07,
    n_iterations = 100
)

"hvf" column was found in the feats metadata information and will be
 used to select highly variable features

Setting dimension reduction [cell][rna] pca

Setting dimension reduction [cell][rna] tsne

Setting dimension reduction [cell][rna] umap

options(repr.plot.width = 11, repr.plot.height = 10)

cosmx@cell_metadata$cell$rna$leiden_clus <- as.character(cosmx@cell_metadata$cell$rna$leiden_clus)
SpatDimPlot(
    cosmx,
    image = TRUE,  # use all fovs
    group_by = "leiden_clus",
    points_size = 1
)

# Small Subset Visualization

options(repr.plot.width = 13, repr.plot.height = 6)

#subset a Giotto object based on spatial locations
smallfov <- subsetGiottoLocs(cosmx,
    x_max = 3000,
    x_min = 1000,
    y_max = -157800,
    y_min = -159800
)

#extract all genes observed in new object
smallfeats <- fDataDT(smallfov)[, feat_ID]

#plot all genes
ext <- c(900, 3100, -159900, -157700)
p1 <- SpatDimPlot(
    smallfov,
    image = "composite_fov002",
    ext = ext,
    shapes = TRUE,
    features = smallfeats,
    points_size = 0.1,
    shapes_border_size = 0.1,
    shapes_border_color = "white",
    legend.position = "none"
)
p2 <- SpatDimPlot(
    smallfov,
    image = "composite_fov002",
    ext = ext,
    shapes = TRUE,
    points_size = 0.1,
    shapes_border_size = 0.1,
    shapes_border_color = "white",
    shapes_alpha = 0
)
p1 + p2

options(repr.plot.width = 7, repr.plot.height = 6)

SpatDimPlot(
    smallfov,
    image = "composite_fov002",
    ext = ext,
    shapes = TRUE,
    points_size = 0.1,
    shapes_fill_by = "leiden_clus",
    shapes_border_size = 0.1,
    shapes_border_color = "white"
)

# create spatial network based on physical distance of cell centroids
cosmx <- createSpatialNetwork(cosmx,
    minimum_k = 2,
    maximum_distance_delaunay = 50
)

# perform Binary Spatial Extraction of genes - NOTE: Depending on your system this could take time
km_spatialfeats <- binSpect(cosmx)

Setting spatial network [cell] Delaunay_network


This is the single parameter version of binSpect


1. matrix binarization complete


2. spatial enrichment test completed


3. (optional) average expression of high
 expressing cells calculated


4. (optional) number of high expressing cells
 calculated

options(repr.plot.width = 14, repr.plot.height = 20)

SpatFeaturePlot(
    cosmx,
    features = km_spatialfeats$feats[1:10],
    points_border_size = 0.01,
    points_size = 0.01,
    facet_ncol = 3,
    padding = 0
)

x <- sessionInfo()
x <- capture.output(print(x))
# hide the BLAS/LAPACK paths
x <- x[!startsWith(x, "BLAS/LAPACK:")]
cat(paste(x, collapse = "\n"))

R version 4.4.3 (2025-02-28)
Platform: x86_64-conda-linux-gnu
Running under: Red Hat Enterprise Linux 8.10 (Ootpa)

Matrix products: default

locale:
 [1] LC_CTYPE=en_US.UTF-8       LC_NUMERIC=C              
 [3] LC_TIME=en_US.UTF-8        LC_COLLATE=C              
 [5] LC_MONETARY=en_US.UTF-8    LC_MESSAGES=en_US.UTF-8   
 [7] LC_PAPER=en_US.UTF-8       LC_NAME=C                 
 [9] LC_ADDRESS=C               LC_TELEPHONE=C            
[11] LC_MEASUREMENT=en_US.UTF-8 LC_IDENTIFICATION=C       

time zone: America/Chicago
tzcode source: system (glibc)

attached base packages:
[1] stats     graphics  grDevices utils     datasets  methods   base     

other attached packages:
[1] scplotter_0.4.0   Giotto_4.2.1      GiottoClass_0.4.7

loaded via a namespace (and not attached):
  [1] fs_1.6.6                    matrixStats_1.5.0          
  [3] spatstat.sparse_3.1-0       bitops_1.0-9               
  [5] devtools_2.4.5              httr_1.4.7                 
  [7] RColorBrewer_1.1-3          repr_1.1.7                 
  [9] profvis_0.4.0               tools_4.4.3                
 [11] sctransform_0.4.2           backports_1.5.0            
 [13] R6_2.6.1                    uwot_0.2.3                 
 [15] lazyeval_0.2.2              urlchecker_1.0.1           
 [17] withr_3.0.2                 sp_2.2-0                   
 [19] gridExtra_2.3               GiottoUtils_0.2.5          
 [21] progressr_0.15.1            quantreg_6.00              
 [23] cli_3.6.5                   Biobase_2.62.0             
 [25] spatstat.explore_3.4-3      fastDummies_1.7.5          
 [27] iNEXT_3.0.1                 Seurat_5.3.0               
 [29] spatstat.data_3.1-6         ggridges_0.5.6             
 [31] pbapply_1.7-2               pbdZMQ_0.3-14              
 [33] stringdist_0.9.15           parallelly_1.45.0          
 [35] sessioninfo_1.2.3           VGAM_1.1-13                
 [37] rstudioapi_0.17.1           generics_0.1.4             
 [39] shape_1.4.6.1               gtools_3.9.5               
 [41] ica_1.0-3                   spatstat.random_3.4-1      
 [43] dplyr_1.1.4                 Matrix_1.7-3               
 [45] S4Vectors_0.40.2            abind_1.4-5                
 [47] terra_1.8-42                lifecycle_1.0.4            
 [49] SummarizedExperiment_1.32.0 SparseArray_1.2.4          
 [51] Rtsne_0.17                  grid_4.4.3                 
 [53] promises_1.3.2              crayon_1.5.3               
 [55] miniUI_0.1.2                lattice_0.22-7             
 [57] cowplot_1.1.3               pillar_1.10.2              
 [59] GenomicRanges_1.54.1        rjson_0.2.23               
 [61] future.apply_1.20.0         codetools_0.2-20           
 [63] glue_1.8.0                  spatstat.univar_3.1-3      
 [65] data.table_1.17.4           remotes_2.5.0              
 [67] vctrs_0.6.5                 png_0.1-8                  
 [69] spam_2.11-1                 gtable_0.3.6               
 [71] assertthat_0.2.1            cachem_1.1.0               
 [73] S4Arrays_1.2.1              mime_0.13                  
 [75] tidygraph_1.3.0             survival_3.8-3             
 [77] SingleCellExperiment_1.24.0 ellipsis_0.3.2             
 [79] scRepertoire_2.2.1          fitdistrplus_1.2-2         
 [81] ROCR_1.0-11                 nlme_3.1-168               
 [83] usethis_3.1.0               RcppAnnoy_0.0.22           
 [85] evd_2.3-7.1                 GenomeInfoDb_1.38.8        
 [87] rprojroot_2.0.4             irlba_2.3.5.1              
 [89] KernSmooth_2.23-26          plotthis_0.7.0             
 [91] colorspace_2.1-1            BiocGenerics_0.48.1        
 [93] tidyselect_1.2.1            compiler_4.4.3             
 [95] SparseM_1.84-2              xml2_1.3.8                 
 [97] desc_1.4.3                  ggdendro_0.2.0             
 [99] DelayedArray_0.28.0         plotly_4.10.4              
[101] checkmate_2.3.2             scales_1.4.0               
[103] lmtest_0.9-40               rappdirs_0.3.3             
[105] stringr_1.5.1               digest_0.6.37              
[107] goftest_1.2-3               spatstat.utils_3.1-4       
[109] XVector_0.42.0              htmltools_0.5.8.1          
[111] GiottoVisuals_0.2.12        pkgconfig_2.0.3            
[113] base64enc_0.1-3             MatrixGenerics_1.14.0      
[115] fastmap_1.2.0               rlang_1.1.6                
[117] GlobalOptions_0.1.2         htmlwidgets_1.6.4          
[119] shiny_1.10.0                farver_2.1.2               
[121] zoo_1.8-14                  jsonlite_2.0.0             
[123] RCurl_1.98-1.17             magrittr_2.0.3             
[125] GenomeInfoDbData_1.2.11     dotCall64_1.2              
[127] patchwork_1.3.0             IRkernel_1.3.2             
[129] Rcpp_1.0.14                 evmix_2.12                 
[131] ggnewscale_0.5.1            viridis_0.6.5              
[133] reticulate_1.42.0           truncdist_1.0-2            
[135] stringi_1.8.7               ggalluvial_0.12.5          
[137] ggraph_2.2.1                zlibbioc_1.48.2            
[139] MASS_7.3-64                 plyr_1.8.9                 
[141] pkgbuild_1.4.8              parallel_4.4.3             
[143] listenv_0.9.1               ggrepel_0.9.6              
[145] forcats_1.0.0               deldir_2.0-4               
[147] graphlayouts_1.2.2          IRdisplay_1.1              
[149] splines_4.4.3               gridtext_0.1.5             
[151] tensor_1.5                  circlize_0.4.16            
[153] colorRamp2_0.1.0            igraph_2.0.3               
[155] uuid_1.2-1                  spatstat.geom_3.4-1        
[157] cubature_2.1.4              RcppHNSW_0.6.0             
[159] reshape2_1.4.4              stats4_4.4.3               
[161] pkgload_1.4.0               evaluate_1.0.3             
[163] SeuratObject_5.1.0          tweenr_2.0.3               
[165] httpuv_1.6.15               MatrixModels_0.5-4         
[167] RANN_2.6.2                  tidyr_1.3.1                
[169] purrr_1.0.4                 polyclip_1.10-7            
[171] future_1.58.0               scattermore_1.2            
[173] ggplot2_3.5.2               ggforce_0.4.2              
[175] xtable_1.8-4                RSpectra_0.16-2            
[177] later_1.4.2                 viridisLite_0.4.2          
[179] gsl_2.1-8                   tibble_3.2.1               
[181] memoise_2.0.1               IRanges_2.36.0             
[183] cluster_2.1.8.1             globals_0.18.0

scplotter to work with Nanostring CosMx data prepared by Giotto¶