diff --git a/DESCRIPTION b/DESCRIPTION
index 8cff1d619..5bd133691 100644
--- a/DESCRIPTION
+++ b/DESCRIPTION
@@ -1,6 +1,6 @@
 Package: Giotto
 Title: Spatial Single-Cell Transcriptomics Toolbox
-Version: 4.1.1
+Version: 4.1.2
 Authors@R: c(
 	person("Ruben", "Dries", email = "rubendries@gmail.com",
 		 role = c("aut", "cre"), comment = c(ORCID = "0000-0001-7650-7754")),
@@ -36,18 +36,14 @@ Imports:
     BiocParallel,
     BiocSingular,
     checkmate,
-    cowplot (>= 0.9.4),
     data.table (>= 1.12.2),
     dbscan (>= 1.1-3),
     ggplot2 (>= 3.1.1),
-    GiottoUtils (>= 0.1.9),
-    GiottoVisuals (>= 0.2.4),
+    GiottoUtils (>= 0.1.12),
+    GiottoVisuals (>= 0.2.5),
     igraph (>= 1.2.4.1),
-    jsonlite,
     Matrix (>= 1.6-2),
     MatrixGenerics,
-    progressr,
-    reshape2,
     reticulate (>= 1.25),
     scales (>= 1.0.0),
     sparseMatrixStats,
diff --git a/NAMESPACE b/NAMESPACE
index f9eddf16a..559afbf2c 100644
--- a/NAMESPACE
+++ b/NAMESPACE
@@ -14,7 +14,6 @@ export("prov<-")
 export("sankeyLabel<-")
 export("sankeyRelate<-")
 export("spatUnit<-")
-export(PAGEEnrich)
 export(activeFeatType)
 export(activeSpatUnit)
 export(addCellIntMetadata)
@@ -138,7 +137,6 @@ export(createSpatLocsObj)
 export(createSpatNetObj)
 export(createSpatialDefaultGrid)
 export(createSpatialDelaunayNetwork)
-export(createSpatialEnrich)
 export(createSpatialFeaturesKNNnetwork)
 export(createSpatialGenomicsObject)
 export(createSpatialGrid)
@@ -262,10 +260,8 @@ export(giottoToSeuratV4)
 export(giottoToSeuratV5)
 export(giottoToSpatialExperiment)
 export(heatmSpatialCorFeats)
-export(heatmSpatialCorGenes)
 export(hexVertices)
 export(hist)
-export(hyperGeometricEnrich)
 export(identifyTMAcores)
 export(importCosMx)
 export(importVisiumHD)
@@ -339,7 +335,6 @@ export(print.combIcfObject)
 export(print.icfObject)
 export(processGiotto)
 export(prov)
-export(rankEnrich)
 export(rankSpatialCorGroups)
 export(read10xAffineImage)
 export(readCellMetadata)
@@ -375,7 +370,6 @@ export(runGiottoHarmony)
 export(runHyperGeometricEnrich)
 export(runIntegratedUMAP)
 export(runPAGEEnrich)
-export(runPAGEEnrich_OLD)
 export(runPCA)
 export(runPCAprojection)
 export(runPCAprojectionBatch)
diff --git a/NEWS.md b/NEWS.md
index 5269365d9..b62e77cec 100644
--- a/NEWS.md
+++ b/NEWS.md
@@ -1,9 +1,27 @@
 
+# Giotto 4.1.2
+
+## Breaking changes
+* remove deprecated `PAGEEnrich()`. Use `runPAGEEnrich()` instead
+* remove deprecated `rankEnrich()`. Use `runRankEnrich()` instead
+* remove deprecated `hyperGeometricEnrich()`. Use `runHyperGeometricEnrich()` instead
+* remove deprecated `createSpatialEnrich()`. Use `runSpatialEnrich()` instead
+* remove deprecated `heatmSpatialCorGenes()`. Use `heatmSpatialCorFeats()` instead
+* remove deprecated `runPAGEEnrich_OLD()`. Use `runPAGEEnrich()` instead
+
+## Changes
+* `limma`, `plotly`, and `Rtsne` moved to Suggests
+* move `progressr` and `jsonlite` dependencies to GiottoUtils v0.1.12
+* remove `reshape2` dependency.
+
 # Giotto 4.1.1
 
 ## Bug fixes
 * Allow `giottoInstructions` passing for Xenium convenience functions
 
+## Changes
+* Deprecate `screePlot()` `name` in favor of `dim_reduction_name` param
+
 
 # Giotto 4.1.0 (2024/07/31)
 
diff --git a/R/ONTraC_wrapper.R b/R/ONTraC_wrapper.R
index e33abeac7..bf11d6075 100644
--- a/R/ONTraC_wrapper.R
+++ b/R/ONTraC_wrapper.R
@@ -398,7 +398,7 @@ plotNicheClusterConnectivity <- function( # nolint: object_name_linter.
   )
   gpl <- gpl + ggraph::scale_edge_alpha(range = c(0.1, 1))
   gpl <- gpl + ggraph::scale_edge_colour_gradientn(
-    colours = RColorBrewer::brewer.pal(9, "Reds"),
+    colours = getColors("Reds", 9, src = "RColorBrewer"),
     name = "Value"
   )
 
@@ -764,4 +764,4 @@ plotCellTypeNTScore <- function(gobject, # nolint: object_name_linter.
     save_param = save_param,
     else_return = NULL
   ))
-}
\ No newline at end of file
+}
diff --git a/R/auxiliary_giotto.R b/R/auxiliary_giotto.R
index 3276c90f5..8fbb26857 100644
--- a/R/auxiliary_giotto.R
+++ b/R/auxiliary_giotto.R
@@ -933,19 +933,12 @@ filterGiotto <- function(
     )
 
     ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ###
-    gobject <- set_expression_values(
-        gobject = gobject,
-        values = norm_expr
-    )
-
-    gobject <- set_expression_values(
-        gobject = gobject,
-        values = norm_scaled_expr
-    )
+    gobject <- setGiotto(gobject, norm_expr, initialize = FALSE)
+    gobject <- setGiotto(gobject, norm_scaled_expr, initialize = FALSE)
     ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ###
 
     ## 6. return Giotto object
-    return(gobject)
+    return(initialize(gobject))
 }
 
 
@@ -991,13 +984,9 @@ filterGiotto <- function(
     )
 
     ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ###
-    gobject <- set_expression_values(
-        gobject = gobject,
-        values = norm_feats_cells
-    )
+    gobject <- setGiotto(giotto, norm_feats_cells)
     ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ###
 
-
     return(gobject)
 }
 
@@ -1088,10 +1077,7 @@ filterGiotto <- function(
     )
 
     ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ###
-    gobject <- set_expression_values(
-        gobject = gobject,
-        values = z
-    )
+    gobject <- setGiotto(gobject, z)
     ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ###
 
     return(gobject)
@@ -1391,10 +1377,7 @@ adjustGiottoMatrix <- function(
         )
 
         ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ###
-        gobject <- set_expression_values(
-            gobject = gobject,
-            values = adjusted_matrix
-        )
+        gobject <- setGiotto(gobject, adjusted_matrix)
         ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ###
 
         ## update parameters used ##
@@ -1456,14 +1439,14 @@ processGiotto <- function(
         adjust_params = list(),
         verbose = TRUE) {
     # filter Giotto
-    if (verbose == TRUE) message("1. start filter step")
+    vmsg(.v = verbose, "1. start filter step")
     if (!inherits(filter_params, "list")) {
         stop("filter_params need to be a list of parameters for filterGiotto")
     }
     gobject <- do.call("filterGiotto", c(gobject = gobject, filter_params))
 
     # normalize Giotto
-    if (verbose == TRUE) message("2. start normalization step")
+    vmsg(.v = verbose, "2. start normalization step")
     if (!inherits(norm_params, "list")) {
         stop("norm_params need to be a list of parameters for normalizeGiotto")
     }
@@ -1623,10 +1606,7 @@ addFeatStatistics <- function(
                 "mean_expr_det"
             ) := NULL]
             ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ###
-            gobject <- set_feature_metadata(gobject,
-                metadata = feat_metadata,
-                verbose = FALSE
-            )
+            gobject <- setGiotto(gobject, feat_metadata, verbose = FALSE)
             ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ###
         }
 
diff --git a/R/clustering.R b/R/clustering.R
index 65c36846c..c0f767dbf 100644
--- a/R/clustering.R
+++ b/R/clustering.R
@@ -1223,8 +1223,8 @@ doSNNCluster <- function(
 #' @param set_seed set seed (default = TRUE)
 #' @param seed_number number for seed
 #' @returns if return_gobject = TRUE: giotto object with new clusters appended to cell metadata
-#' @details The default settings will use dimension reduction results as input. 
-#' Set dim_reduction_to_use = NULL if you want to directly use expression values as input. 
+#' @details The default settings will use dimension reduction results as input.
+#' Set dim_reduction_to_use = NULL if you want to directly use expression values as input.
 #' By providing a feature vector to feats_to_use you can subset the expression matrix.
 #' @seealso  \code{\link[stats]{kmeans}}
 #' @examples
@@ -1254,7 +1254,7 @@ doKmeans <- function(
         return_gobject = TRUE,
         set_seed = TRUE,
         seed_number = 1234) {
-  
+
     # Set feat_type and spat_unit
     spat_unit <- set_default_spat_unit(
         gobject = gobject,
@@ -1276,9 +1276,9 @@ doKmeans <- function(
 
     ## using dimension reduction ##
     if(!is.null(dim_reduction_to_use)) {
-      
+
       # use only available dimensions if dimensions < dimensions_to_use
-      dim_coord <- get_dimReduction(
+      dim_coord <- getDimReduction(
         gobject = gobject,
         spat_unit = spat_unit,
         feat_type = feat_type,
@@ -1287,14 +1287,14 @@ doKmeans <- function(
         name = dim_reduction_name,
         output = "dimObj"
       )
-      
+
       dimensions_to_use <- dimensions_to_use[
         dimensions_to_use %in% seq_len(ncol(dim_coord[]))
       ]
       matrix_to_use <- dim_coord[][, dimensions_to_use]
-      
+
     } else {
-      
+
       ## using original matrix ##
       expr_values <- getExpression(
         gobject = gobject,
@@ -1303,21 +1303,21 @@ doKmeans <- function(
         values = expression_values,
         output = "exprObj"
       )
-      
+
       # subset expression matrix
       if (!is.null(feats_to_use)) {
         expr_values[] <- expr_values[][
           rownames(expr_values[]) %in% feats_to_use,
         ]
       }
-      
+
       # features as columns
       # cells as rows
       matrix_to_use <- t_flex(expr_values[])
-      
+
     }
-    
-    
+
+
     ## distance
     if (distance_method == "original") {
         celldist <- matrix_to_use
@@ -2817,7 +2817,7 @@ getClusterSimilarity <- function(
 
     # correlation matrix
     cormatrix <- cor_flex(x = testmatrix, method = cor)
-    cor_table <- data.table::as.data.table(reshape2::melt(cormatrix))
+    cor_table <- melt_matrix(cormatrix)
     data.table::setnames(
         cor_table,
         old = c("Var1", "Var2"), c("group1", "group2")
diff --git a/R/convenience_cosmx.R b/R/convenience_cosmx.R
index 4a107002f..c3a503a6a 100644
--- a/R/convenience_cosmx.R
+++ b/R/convenience_cosmx.R
@@ -723,8 +723,8 @@ setMethod("$<-", signature("CosmxReader"), function(x, name, value) {
     )
 
     fovs <- fovs %null% .cosmx_imgname_fovparser(path) # ALL if NULL
-    progressr::with_progress({
-        p <- progressr::progressor(along = fovs)
+    with_pbar({
+        p <- pbar(along = fovs)
 
         gpolys <- lapply(fovs, function(f) {
             segfile <- Sys.glob(paths = sprintf("%s/*F%03d*", path, f))
@@ -925,8 +925,8 @@ setMethod("$<-", signature("CosmxReader"), function(x, name, value) {
     fovs <- fovs %null% .cosmx_imgname_fovparser(path) # ALL if NULL
     verbose <- verbose %null% TRUE
 
-    progressr::with_progress({
-        p <- progressr::progressor(along = fovs)
+    with_pbar({
+        p <- pbar(along = fovs)
 
         gimg_list <- lapply(fovs, function(f) {
             imgfile <- Sys.glob(paths = sprintf("%s/*F%03d*", path, f))
@@ -1645,7 +1645,7 @@ createGiottoCosMxObject <- function(
         provenance = "cell_agg"
     )
 
-    cosmx_gobject <- set_expression_values(cosmx_gobject, values = s4_expr)
+    cosmx_gobject <- setGiotto(cosmx_gobject, s4_expr)
 
     # Add spatial locations
     if (isTRUE(verbose)) wrap_msg(
@@ -1669,9 +1669,8 @@ createGiottoCosMxObject <- function(
         provenance = "cell_agg"
     )
 
-    cosmx_gobject <- set_spatial_locations(cosmx_gobject, spatlocs = locsObj)
-    cosmx_gobject <- set_spatial_locations(cosmx_gobject,
-                                           spatlocs = locsObj_fov)
+    cosmx_gobject <- setGiotto(cosmx_gobject, locsObj, initialize = FALSE)
+    cosmx_gobject <- setGiotto(cosmx_gobject, locsObj_fov, initialize = FALSE)
 
     # initialize cell and feat IDs and metadata slots for 'cell_agg' spat_unit
     agg_cell_ID <- colnames(s4_expr[])
diff --git a/R/convenience_general.R b/R/convenience_general.R
index 1d626c5bf..d5815eaa4 100644
--- a/R/convenience_general.R
+++ b/R/convenience_general.R
@@ -656,7 +656,7 @@ addVisiumPolygons <- function(gobject,
         return(NULL)
     }
 
-    json_scalefactors <- jsonlite::read_json(json_path)
+    json_scalefactors <- read_json(json_path)
 
     # Intial assertion that json dimensions are appropriate
     checkmate::assert_list(
diff --git a/R/convenience_visiumHD.R b/R/convenience_visiumHD.R
index a0c3ca5fe..645b76ee4 100644
--- a/R/convenience_visiumHD.R
+++ b/R/convenience_visiumHD.R
@@ -746,7 +746,7 @@ setMethod("$<-", signature("VisiumHDReader"), function(x, name, value) {
     json_path = file.path(path, 'scalefactors_json.json')
     checkmate::assert_file_exists(json_path)
 
-    json_scalefactors <- jsonlite::read_json(json_path)
+    json_scalefactors <- read_json(json_path)
 
     expected_json_names <- c(
         "regist_target_img_scalef", # NEW as of 2023
diff --git a/R/convenience_xenium.R b/R/convenience_xenium.R
index e8b437286..6fc5842b7 100644
--- a/R/convenience_xenium.R
+++ b/R/convenience_xenium.R
@@ -193,7 +193,7 @@ setMethod(
         # decide micron scaling
         if (length(obj@micron) == 0) { # if no value already set
             if (!is.null(experiment_info_path)) {
-                obj@micron <- jsonlite::fromJSON(
+                obj@micron <- fromJSON(
                     experiment_info_path)$pixel_size
             } else {
                 warning(wrap_txt("No .xenium file found.
@@ -914,7 +914,7 @@ importXenium <- function(
     gene_ids <- match.arg(gene_ids, c("symbols", "ensembl"))
 
     # tested on v1.6
-    j <- jsonlite::fromJSON(path)
+    j <- fromJSON(path)
     # j$metadata # dataset meta
     # j$payload # main content
     # j$payload$chemistry # panel chemistry used
@@ -1102,8 +1102,8 @@ importXenium <- function(
     # micron
     checkmate::assert_numeric(micron)
 
-    progressr::with_progress({
-        p <- progressr::progressor(along = path)
+    with_pbar({
+        p <- pbar(along = path)
 
         gimg_list <- lapply(seq_along(path), function(img_i) {
             gimg <- .xenium_image_single(
diff --git a/R/differential_expression.R b/R/differential_expression.R
index 5eead79a0..d84d06f53 100644
--- a/R/differential_expression.R
+++ b/R/differential_expression.R
@@ -279,8 +279,8 @@ findScranMarkers_one_vs_all <- function(
 
 
     # save list
-    progressr::with_progress({
-        pb <- progressr::progressor(along = uniq_clusters)
+    with_pbar({
+        pb <- pbar(along = uniq_clusters)
         result_list <- lapply(
             seq_along(uniq_clusters),
             function(clus_i) {
@@ -734,8 +734,8 @@ findGiniMarkers_one_vs_all <- function(
 
 
     # GINI
-    progressr::with_progress({
-        pb <- progressr::progressor(along = uniq_clusters)
+    with_pbar({
+        pb <- pbar(along = uniq_clusters)
         result_list <- lapply(
             seq_along(uniq_clusters),
             function(clus_i) {
diff --git a/R/dimension_reduction.R b/R/dimension_reduction.R
index 7b4f8bc18..de7906d1c 100644
--- a/R/dimension_reduction.R
+++ b/R/dimension_reduction.R
@@ -541,9 +541,7 @@ runPCA <- function(
         )
 
         ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ###
-        gobject <- set_dimReduction(
-            gobject = gobject, dimObject = dimObject, verbose = verbose
-        )
+        gobject <- setGiotto(gobject, dimObject, verbose = verbose)
         ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ###
 
 
@@ -945,7 +943,7 @@ runPCAprojection <- function(
         )
 
         ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ###
-        gobject <- set_dimReduction(gobject = gobject, dimObject = dimObject)
+        gobject <- setGiotto(gobject, dimObject, verbose = verbose)
         ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ###
 
 
@@ -1349,7 +1347,7 @@ runPCAprojectionBatch <- function(
         )
 
         ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ###
-        gobject <- set_dimReduction(gobject = gobject, dimObject = dimObject)
+        gobject <- setGiotto(gobject, dimObject, verbose = verbose)
         ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ###
 
 
@@ -1376,7 +1374,8 @@ runPCAprojectionBatch <- function(
 #' @inheritParams data_access_params
 #' @inheritParams plot_output_params
 #' @inheritParams create_screeplot
-#' @param name name of PCA object if available
+#' @param dim_reduction_name name of PCA
+#' @param name deprecated
 #' @param expression_values expression values to use
 #' @param reduction cells or features
 #' @param method which implementation to use
@@ -1402,7 +1401,8 @@ screePlot <- function(
         gobject,
         spat_unit = NULL,
         feat_type = NULL,
-        name = NULL,
+        dim_reduction_name = NULL,
+        name = deprecated(),
         expression_values = c("normalized", "scaled", "custom"),
         reduction = c("cells", "feats"),
         method = c("irlba", "exact", "random", "factominer"),
@@ -1419,6 +1419,17 @@ screePlot <- function(
         save_param = list(),
         default_save_name = "screePlot",
         ...) {
+
+    if (is_present(name)) {
+        deprecate_warn(
+            when = "4.1.1",
+            what = "screePlot(name = )",
+            with = "screePlot(dim_reduction_name = )"
+        )
+    } else {
+        name <- dim_reduction_name # shorter varname
+    }
+
     # Set feat_type and spat_unit
     spat_unit <- set_default_spat_unit(
         gobject = gobject,
@@ -1629,7 +1640,7 @@ create_screeplot <- function(eigs, ncp = 20, ylim = c(0, 20)) {
     savelist <- list(pl, cpl)
 
     ## combine plots with cowplot
-    combo_plot <- cowplot::plot_grid(
+    combo_plot <- plot_grid(
         plotlist = savelist,
         ncol = 1,
         rel_heights = c(1),
@@ -2155,7 +2166,7 @@ runUMAP <- function(
         ## using dimension reduction ##
         if (!is.null(dim_reduction_to_use)) {
             ## TODO: check if reduction exists
-            dimObj_to_use <- get_dimReduction(
+            dimObj_to_use <- getDimReduction(
                 gobject = gobject,
                 spat_unit = spat_unit,
                 feat_type = feat_type,
@@ -2278,10 +2289,7 @@ runUMAP <- function(
 
 
             ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ###
-            gobject <- set_dimReduction(
-                gobject = gobject,
-                dimObject = dimObject
-            )
+            gobject <- setGiotto(gobject, dimObject, verbose = verbose)
             ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ###
 
 
@@ -2553,7 +2561,7 @@ runUMAPprojection <- function(
 
 
         ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ###
-        gobject <- set_dimReduction(gobject = gobject, dimObject = dimObject)
+        gobject <- setGiotto(gobject, dimObject, verbose = verbose)
         ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ###
 
         ## update parameters used ##
@@ -2781,10 +2789,7 @@ runtSNE <- function(
             )
 
             ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ###
-            gobject <- set_dimReduction(
-                gobject = gobject,
-                dimObject = dimObject
-            )
+            gobject <- setGiotto(gobject, dimObject, verbose = verbose)
             ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ###
 
             ## update parameters used ##
@@ -2918,7 +2923,7 @@ runGiottoHarmony <- function(
     ## using dimension reduction ##
     if (!is.null(dim_reduction_to_use)) {
         ## TODO: check if reduction exists
-        matrix_to_use <- get_dimReduction(
+        matrix_to_use <- getDimReduction(
             gobject = gobject,
             spat_unit = spat_unit,
             feat_type = feat_type,
@@ -2987,14 +2992,14 @@ runGiottoHarmony <- function(
     colnames(harmony_results) <- paste0("Dim.", seq_len(ncol(harmony_results)))
     rownames(harmony_results) <- rownames(matrix_to_use)
 
-    harmdimObject <- create_dim_obj(
+    harmdimObject <- createDimObj(
+        coordinates = harmony_results,
         name = name,
         spat_unit = spat_unit,
         feat_type = feat_type,
-        provenance = provenance,
+        method = "harmony",
         reduction = "cells", # set to spat_unit?
-        reduction_method = "harmony",
-        coordinates = harmony_results,
+        provenance = provenance,
         misc = NULL
     )
 
@@ -3016,10 +3021,7 @@ runGiottoHarmony <- function(
         }
 
         ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ###
-        gobject <- set_dimReduction(
-            gobject = gobject,
-            dimObject = harmdimObject
-        )
+        gobject <- setGiotto(gobject, harmdimObject, verbose = verbose)
         ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ###
 
 
diff --git a/R/general_help.R b/R/general_help.R
index c10ae4ffb..cc0b3daed 100644
--- a/R/general_help.R
+++ b/R/general_help.R
@@ -744,7 +744,7 @@ read10xAffineImage <- function(
     checkmate::assert_file_exists(imagealignment_path)
     if (!is.numeric(micron)) {
         checkmate::assert_file_exists(micron)
-        micron <- jsonlite::read_json(micron)$pixel_size
+        micron <- read_json(micron)$pixel_size
     }
 
     aff <- data.table::fread(imagealignment_path) %>%
@@ -859,8 +859,8 @@ readPolygonFilesVizgenHDF5_old <- function(
 
     # append data from all FOVs to single list
     init <- proc.time()
-    progressr::with_progress({
-        pb <- progressr::progressor(along = hdf5_boundary_selected_list)
+    with_pbar({
+        pb <- pbar(along = hdf5_boundary_selected_list)
         read_list <- lapply_flex(seq_along(hdf5_boundary_selected_list),
             cores = cores,
             future.packages = c("rhdf5", "Rhdf5lib"),
@@ -933,8 +933,8 @@ readPolygonFilesVizgenHDF5_old <- function(
 
 
     # create Giotto polygons and add them to gobject
-    progressr::with_progress({
-        pb <- progressr::progressor(along = result_list_rbind)
+    with_pbar({
+        pb <- pbar(along = result_list_rbind)
         smooth_cell_polygons_list <- lapply_flex(seq_along(result_list_rbind),
             cores = cores, function(i) {
                 dfr_subset <- result_list_rbind[[i]][, .(x, y, cell_id)]
@@ -1074,8 +1074,8 @@ readPolygonFilesVizgenHDF5 <- function(
 
     # append data from all FOVs to single list
     init <- Sys.time()
-    progressr::with_progress({
-        pb <- progressr::progressor(length(hdf5_boundary_selected_list) / 5)
+    with_pbar({
+        pb <- pbar(length(hdf5_boundary_selected_list) / 5)
         read_list <- lapply_flex(seq_along(hdf5_boundary_selected_list),
             future.packages = c("rhdf5", "Rhdf5lib"),
             function(init, z_indices, segm_to_use, bound_i) {
@@ -1167,8 +1167,8 @@ readPolygonFilesVizgenHDF5 <- function(
 
     # **** sequential method ****
     if (!isTRUE(create_gpoly_parallel)) {
-        progressr::with_progress({
-            pb <- progressr::progressor(along = z_read_DT)
+        with_pbar({
+            pb <- pbar(along = z_read_DT)
             smooth_cell_polygons_list <- lapply(
                 seq_along(z_read_DT), function(i) {
                     dfr_subset <- z_read_DT[[i]][, .(x, y, cell_id)]
@@ -1213,8 +1213,8 @@ readPolygonFilesVizgenHDF5 <- function(
     # **** parallel methods ****
     # no binning
     if (!is.numeric(create_gpoly_bin)) {
-        progressr::with_progress({
-            pb <- progressr::progressor(along = z_read_DT)
+        with_pbar({
+            pb <- pbar(along = z_read_DT)
             smooth_cell_polygons_list <- lapply_flex(
                 seq_along(z_read_DT),
                 future.packages = c("terra", "stats", "data.table"),
@@ -1282,8 +1282,8 @@ readPolygonFilesVizgenHDF5 <- function(
 
         bin_steps <- sum(unlist(lapply(dfr_subset, length)))
 
-        progressr::with_progress({
-            pb <- progressr::progressor(steps = bin_steps)
+        with_pbar({
+            pb <- pbar(steps = bin_steps)
             smooth_cell_polygons_list <- lapply( # sequential across z index
                 seq_along(dfr_subset),
                 function(i) {
diff --git a/R/giotto_viewer.R b/R/giotto_viewer.R
index c21f66f75..8b044855d 100644
--- a/R/giotto_viewer.R
+++ b/R/giotto_viewer.R
@@ -212,9 +212,10 @@ exportGiottoViewer <- function(
 
     # data.table variables
     sdimx <- sdimy <- NULL
-    spatial_location <- get_spatial_locations(
+    spatial_location <- getSpatialLocations(
         gobject = gobject,
-        spat_loc_name = spat_loc_name
+        spat_unit = spat_unit,
+        name = spat_loc_name
     )
     spatial_location <- spatial_location[, .(sdimx, sdimy)]
     write.table(spatial_location,
@@ -241,6 +242,7 @@ exportGiottoViewer <- function(
         cell_metadata <- combineMetadata(
             gobject = gobject,
             feat_type = feat,
+            spat_unit = spat_unit,
             spat_enr_names = spat_enr_names
         )
 
@@ -361,11 +363,12 @@ exportGiottoViewer <- function(
     values <- match.arg(expression_values, unique(c("scaled", "normalized", "custom", expression_values)))
 
     for (feat in feat_type) {
-        expr_values <- get_expression_values(
+        expr_values <- getExpression(
             gobject = gobject,
             spat_unit = spat_unit,
             feat_type = feat,
-            values = values
+            values = values,
+            output = "matrix"
         )
         expr_values <- as.matrix(expr_values)
 
diff --git a/R/image_registration.R b/R/image_registration.R
index a98e9a280..f163b30f2 100644
--- a/R/image_registration.R
+++ b/R/image_registration.R
@@ -163,7 +163,7 @@
         img_type = img_type
     )))) {
         giottoImage_list <- lapply(
-            X = gobject_list, FUN = get_giottoImage, name = image_unreg,
+            X = gobject_list, FUN = getGiottoImage, name = image_unreg,
             image_type = img_type
         )
         image_corners <- lapply(giottoImage_list, .get_img_corners)
@@ -431,10 +431,10 @@ registerGiottoObjectListFiji <- function(
     spatloc_list <- list()
     for (gobj_i in seq_along(gobject_list)) {
         gobj <- gobject_list[[gobj_i]]
-        spatloc <- get_spatial_locations(
+        spatloc <- getSpatialLocations(
             gobject = gobj,
             spat_unit = spat_unit,
-            spat_loc_name = spatloc_unreg
+            name = spatloc_unreg
         )
         #------ Put all spatial location data together
         spatloc_list[[gobj_i]] <- spatloc
@@ -526,24 +526,24 @@ registerGiottoObjectListFiji <- function(
         # Rename original spatial locations to 'unregistered' if conflicting
         # with output
         if (spatloc_unreg == spatloc_reg_name) {
-            gobj <- set_spatial_locations(
+            gobj <- setSpatialLocations(
                 gobject = gobj,
                 spat_unit = spat_unit,
-                spat_loc_name = spatloc_replace_name,
-                spatlocs = get_spatial_locations(
+                name = spatloc_replace_name,
+                spatlocs = getSpatialLocations(
                     gobject = gobj,
                     spat_unit = spat_unit,
-                    spat_loc_name = spatloc_unreg
+                    name = spatloc_unreg
                 )
             )
         }
 
 
         # Assign registered spatial locations from spatloc_list to gobject_list
-        gobj <- set_spatial_locations(
+        gobj <- setSpatialLocations(
             gobject = gobj,
             spat_unit = spat_unit,
-            spat_loc_name = spatloc_reg_name,
+            name = spatloc_reg_name,
             spatlocs = spatloc_list[[gobj_i]]
         )
 
@@ -637,9 +637,9 @@ registerGiottoObjectListRvision <- function(
     spatloc_list <- list()
     for (gobj_i in seq_along(gobject_list)) {
         gobj <- gobject_list[[gobj_i]]
-        spatloc <- get_spatial_locations(
+        spatloc <- getSpatialLocations(
             gobject = gobj,
-            spat_loc_name = spatloc_unreg,
+            name = spatloc_unreg,
             output = "spatLocsObj",
             copy_obj = TRUE
         )
@@ -744,22 +744,22 @@ registerGiottoObjectListRvision <- function(
         gobj <- gobject_list[[gobj_i]]
         # Rename original spatial locations to 'unregistered'
 
-        unreg_locs <- get_spatial_locations(gobj,
-            spat_loc_name = spatloc_unreg,
+        unreg_locs <- getSpatialLocations(gobj,
+            name = spatloc_unreg,
             copy_obj = FALSE,
             output = "spatLocsObj"
         )
 
         ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ###
-        gobj <- set_spatial_locations(gobj,
+        gobj <- setSpatialLocations(gobj,
             spatlocs = unreg_locs,
-            spat_loc_name = "unregistered"
+            name = "unregistered"
         )
 
         # Assign registered spatial locations from spatloc_list to gobject_list
-        gobj <- set_spatial_locations(gobj,
+        gobj <- setSpatialLocations(gobj,
             spatlocs = spatloc_list[[gobj_i]],
-            spat_loc_name = spatloc_reg_name
+            name = spatloc_reg_name
         )
         ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ###
 
@@ -1033,7 +1033,7 @@ interactiveLandmarkSelection <- function(source, target) {
     GiottoUtils::package_check("shiny")
     GiottoUtils::package_check("ggplot2")
     GiottoUtils::package_check("miniUI")
-    
+
     .create_image_to_plot <- function(x){
         if (inherits(x, "gg")){
             return(x)
@@ -1052,7 +1052,7 @@ interactiveLandmarkSelection <- function(source, target) {
     }
     source_image <- .create_image_to_plot(source)
     target_image <- .create_image_to_plot(target)
-    
+
     # Function to extract the range of x and y values from a ggplot object
     .extract_plot_ranges <- function(plot) {
         data <- ggplot2::ggplot_build(plot)$data[[1]]
@@ -1060,11 +1060,11 @@ interactiveLandmarkSelection <- function(source, target) {
         y_range <- range(data$y, na.rm = TRUE)
         list(x_range = x_range, y_range = y_range)
     }
-    
+
     # Extract ranges for the input plots
     source_ranges <- .extract_plot_ranges(source_image)
     target_ranges <- .extract_plot_ranges(target_image)
-    
+
     ui <- miniUI::miniPage(
         miniUI::gadgetTitleBar("Select Extents and Points"),
         miniUI::miniContentPanel(
@@ -1073,11 +1073,11 @@ interactiveLandmarkSelection <- function(source, target) {
                 shiny::column(6, shiny::plotOutput("plot2", click = "plot2_click"))
             ),
             shiny::fluidRow(
-                shiny::column(6, 
+                shiny::column(6,
                               shiny::sliderInput("xrange1", "X Range for Plot 1", min = source_ranges$x_range[1], max = source_ranges$x_range[2], value = source_ranges$x_range),
                               shiny::sliderInput("yrange1", "Y Range for Plot 1", min = source_ranges$y_range[1], max = source_ranges$y_range[2], value = source_ranges$y_range)
                 ),
-                shiny::column(6, 
+                shiny::column(6,
                               shiny::sliderInput("xrange2", "X Range for Plot 2", min = target_ranges$x_range[1], max = target_ranges$x_range[2], value = target_ranges$x_range),
                               shiny::sliderInput("yrange2", "Y Range for Plot 2", min = target_ranges$y_range[1], max = target_ranges$y_range[2], value = target_ranges$y_range)
                 )
@@ -1092,63 +1092,63 @@ interactiveLandmarkSelection <- function(source, target) {
             )
         )
     )
-    
+
     server <- function(input, output, session) {
         click_history1 <- shiny::reactiveVal(data.frame(x = numeric(), y = numeric()))
         click_history2 <- shiny::reactiveVal(data.frame(x = numeric(), y = numeric()))
-        
+
         output$plot1 <- shiny::renderPlot({
             source_image +
                 ggplot2::coord_cartesian(xlim = input$xrange1, ylim = input$yrange1) +
                 ggplot2::geom_point(data = click_history1(), ggplot2::aes(x = x, y = y), color = "red", size = 4.5)
         })
-        
+
         output$plot2 <- shiny::renderPlot({
             target_image +
                 ggplot2::coord_cartesian(xlim = input$xrange2, ylim = input$yrange2) +
                 ggplot2::geom_point(data = click_history2(), ggplot2::aes(x = x, y = y), color = "blue",size = 4.5)
         })
-        
+
         shiny::observeEvent(input$plot1_click, {
             click <- input$plot1_click
             new_coords <- rbind(click_history1(), data.frame(x = click$x, y = click$y))
             click_history1(new_coords)
         })
-        
+
         shiny::observeEvent(input$plot2_click, {
             click <- input$plot2_click
             new_coords <- rbind(click_history2(), data.frame(x = click$x, y = click$y))
             click_history2(new_coords)
         })
-        
+
         shiny::observeEvent(input$undo1, {
             if (nrow(click_history1()) > 0) {
                 new_coords <- click_history1()[-nrow(click_history1()), , drop = FALSE]
                 click_history1(new_coords)
             }
         })
-        
+
         shiny::observeEvent(input$undo2, {
             if (nrow(click_history2()) > 0) {
                 new_coords <- click_history2()[-nrow(click_history2()), , drop = FALSE]
                 click_history2(new_coords)
             }
         })
-        
+
         output$click_info1 <- shiny::renderPrint({
             click_history1()
         })
-        
+
         output$click_info2 <- shiny::renderPrint({
             click_history2()
         })
-        
+
         shiny::observeEvent(input$done, {
             returnValue <- list(click_history1(),click_history2())
             shiny::stopApp(returnValue)
         })
     }
-    
+
     shiny::runGadget(ui, server)
 }
 
@@ -1184,32 +1184,32 @@ calculateAffineMatrixFromLandmarks <- function(source_df,target_df){
 
 
 #' @name .sift_detect
-#' @title Run SIFT feature detector and descriptor extractor 
+#' @title Run SIFT feature detector and descriptor extractor
 #' @description
 #' Perform feature detector and descriptor extractor on a matrix object or preprocessed image object
 #' @param x input matrix or preprocessed image to extract feature and descriptor from
 #' @param ... additional params to pass to `skimage.feature.SIFT()`
 #' @returns list of keypoints and descriptors
-#' 
+#'
 .sift_detect <- function(x, ..., pkg_ptr) {
-    
+
     if (missing(pkg_ptr)) {
         GiottoUtils::package_check("skimage", repository = "pip:scikit-image")
         SKI <- reticulate::import("skimage", convert = TRUE, delay_load = TRUE)
     } else {
         SKI <- pkg_ptr
     }
-    
+
     # sift object
     SIFT <- SKI$feature$SIFT()
-    
+
     SIFT$detect_and_extract(x)
-    
+
     out <- list(
         keypoints = SIFT$keypoints,
         descriptors = SIFT$descriptors
     )
-    
+
     return(out)
 }
 
@@ -1230,7 +1230,7 @@ calculateAffineMatrixFromLandmarks <- function(source_df,target_df){
 #' Computer Vision, 2004.
 #' @param ... additional params to pass to `skimage.feature.match_descriptors()`
 #' @returns list
-#' 
+#'
 .match_descriptor <- function(
         descriptor_list,
         target_idx = 1L,
@@ -1239,19 +1239,19 @@ calculateAffineMatrixFromLandmarks <- function(source_df,target_df){
         ...,
         pkg_ptr
 ) {
-    
+
     checkmate::assert_list(descriptor_list, min.len = 2L)
     target_idx <- as.integer(target_idx)
-    
+
     if (missing(pkg_ptr)) {
         package_check("skimage", repository = "pip:scikit-image")
         SKI <- reticulate::import("skimage", convert = TRUE, delay_load = TRUE)
     } else {
         SKI <- pkg_ptr
     }
-    
+
     target <- descriptor_list[[target_idx]]
-    
+
     out <- lapply(
         seq_along(descriptor_list),
         function(moving_idx) {
@@ -1263,9 +1263,9 @@ calculateAffineMatrixFromLandmarks <- function(source_df,target_df){
                 ))
                 # directly return all as matches
             }
-            
+
             moving <- descriptor_list[[moving_idx]]
-            
+
             m <- .match_descriptor_single(
                 x = target,
                 y = moving,
@@ -1275,7 +1275,7 @@ calculateAffineMatrixFromLandmarks <- function(source_df,target_df){
             m + 1 # since it is 0 indexed
         }
     )
-    
+
     return(out)
 }
 
@@ -1283,17 +1283,17 @@ calculateAffineMatrixFromLandmarks <- function(source_df,target_df){
 # wrapper for sklearn-image match_descriptors
 # returns a 2 col matrix of x to y index matches
 .match_descriptor_single <- function(x, y,max_ratio, ..., pkg_ptr) {
-    
+
     checkmate::assert_class(x, "matrix")
     checkmate::assert_class(y, "matrix")
-    
+
     if (missing(pkg_ptr)) {
         GiottoUtils::package_check("skimage", repository = "pip:scikit-image")
         SKI <- reticulate::import("skimage", convert = TRUE, delay_load = TRUE)
     } else {
         SKI <- pkg_ptr
     }
-    
+
     match_descriptors <- SKI$feature$match_descriptors
     m <- match_descriptors(
         descriptors1 = x,
@@ -1301,7 +1301,7 @@ calculateAffineMatrixFromLandmarks <- function(source_df,target_df){
         max_ratio = max_ratio,
         ... # max_ratio of 0.6 - 0.8 recommended for sift, cross_check = TRUE
     )
-    
+
     return(m)
 }
 
@@ -1319,7 +1319,7 @@ calculateAffineMatrixFromLandmarks <- function(source_df,target_df){
 #' @param use_single_channel If input is a multichannel image, whether or not to extract single channel, default FALSE
 #' @param single_channel_number Channel number in the multichannel image, required if use_single_channel = TRUE
 #' @returns a matrix array to input to .sift_detect
-#' 
+#'
 #' @export
 preprocessImageToMatrix <- function(x,
                                     invert = F,
@@ -1330,7 +1330,7 @@ preprocessImageToMatrix <- function(x,
                                     use_single_channel = F,
                                     single_channel_number = NULL,
                                     pkg_ptr) {
-    
+
     if (missing(pkg_ptr)) {
         GiottoUtils::package_check("skimage", repository = "pip:scikit-image")
         SKI <- reticulate::import("skimage", convert = TRUE, delay_load = TRUE)
@@ -1339,18 +1339,18 @@ preprocessImageToMatrix <- function(x,
     }
     GiottoUtils::package_check("numpy", repository = "pip:scikit-image")
     np <- reticulate::import("numpy", convert = TRUE, delay_load = TRUE)
-    
+
     image = SKI$io$imread(x)
-    
+
     if (length(dim(image)) >2 & use_single_channel == FALSE){
         image = SKI$color$rgb2gray(image)
-    } 
+    }
     if (use_single_channel  == TRUE) {
         if (is.null(single_channel_number)) {stop("Set use single channel == TRUE, please provide a channel number to continue")}
         image <- image[,,single_channel_number]
     }
-    
-    
+
+
     if (flip_vertical == T){
         image = np$flipud(image)
     }
@@ -1390,20 +1390,20 @@ preprocessImageToMatrix <- function(x,
     } else {
         SKI <- pkg_ptr
     }
-    
+
     # Extract matched keypoints
     src_pts <- keypoints1[match[, 1] + 1, , drop = FALSE]
     dst_pts <- keypoints2[match[, 2] + 1, , drop = FALSE]
-    
+
     estimate_fun <- match.arg(estimate_fun, unique(c('euclidean', 'similarity', 'affine', 'piecewise-affine', 'projective', 'polynomial', estimate_fun)))
-    
-    # Estimate homography matrix 
+
+    # Estimate homography matrix
     ransac_result <- SKI$transform$estimate_transform(
         ttype = estimate_fun,
         src = src_pts,
         dst = dst_pts,
     )
-    
+
     return(ransac_result)
 }
 
@@ -1427,10 +1427,10 @@ preprocessImageToMatrix <- function(x,
     } else {
         SKI <- pkg_ptr
     }
-    
+
     # Ensure the source image array is writable by making a copy
     x_copy <- reticulate::r_to_py(x)$copy()
-    
+
     # Warp the source image to align with the destination image
     warped_image <- SKI$transform$warp(x_copy, model, output_shape = dim(y))
     SKI$io$imsave(outpath,warped_image)
@@ -1455,22 +1455,22 @@ preprocessImageToMatrix <- function(x,
     } else {
         SKI <- pkg_ptr
     }
-    
+
     matplotlib <-reticulate::import("matplotlib", convert = TRUE, delay_load = TRUE)
     np <- reticulate::import("numpy",convert = T, delay_load = T)
     plt <- matplotlib$pyplot
-    
+
     match_py <- reticulate::r_to_py(match)
-    match_py <- np$array(match_py, dtype = np$int32) 
-    
+    match_py <- np$array(match_py, dtype = np$int32)
+
     # Create a subplot
     fig_ax <- plt$subplots(nrows = 1L, ncols = 1L, figsize = c(11, 8))
     fig <- fig_ax[[1]]
     ax <- fig_ax[[2]]
-    
+
     # Plot the matches
     SKI$feature$plot_matches(ax, x, y, keypoints1, keypoints2, match_py, only_matches = TRUE)
-    
+
     ax$axis('off')
     plt$show()
     plt$close()
@@ -1497,31 +1497,31 @@ estimateAutomatedImageRegistrationWithSIFT <- function(x,
                                                        estimate_fun = 'affine',
                                                        save_warp = NULL,
                                                        verbose = T){
-    
+
     GiottoUtils::vmsg(.v = verbose, .is_debug = T,'Detecting features via SIFT... ')
     x_sift <- .sift_detect(x)
     y_sift <- .sift_detect(y)
-    
+
     GiottoUtils::vmsg(.v = verbose, .is_debug = T,'Matching Descriptors via SIFT... ')
     matched <- .match_descriptor_single(x_sift$descriptor, y_sift$descriptor,max_ratio = max_ratio)
-    
+
     if (plot_match == TRUE){
         .plot_matched_descriptors(x, y, x_sift$keypoints, y_sift$keypoints, matched)
     }
-    
-    
+
+
     GiottoUtils::vmsg(.v = verbose, .is_debug = T,'Estimating transformation matrix from matched descriptor... ')
-    estimation <- .estimate_transform_from_matched_descriptor(x_sift$keypoints, 
+    estimation <- .estimate_transform_from_matched_descriptor(x_sift$keypoints,
                                                               y_sift$keypoints,
                                                               matched,
                                                               estimate_fun = estimate_fun)
-    
+
     if (!is.null(save_warp)){
         .warp_transformed_image(x = x,
                                 y = y,
                                 model = estimation$inverse, outpath = save_warp)
     }
-    
+
     return(estimation)
 }
 
diff --git a/R/interactivity.R b/R/interactivity.R
index 27dc3ad5c..2ce831111 100644
--- a/R/interactivity.R
+++ b/R/interactivity.R
@@ -380,7 +380,7 @@ comparePolygonExpression <- function(
     }
 
     # get expression
-    my_expression <- get_expression_values(gobject,
+    my_expression <- getExpression(gobject,
         values = expression_values,
         spat_unit = spat_unit,
         feat_type = feat_type,
@@ -587,7 +587,7 @@ plotPolygons <- function(
     if (is.null(x)) stop("A plot object must be provided")
 
     ## get polygons spatial info
-    polygon_spatVector <- get_polygon_info(
+    polygon_spatVector <- getPolygonInfo(
         gobject = gobject,
         polygon_name = polygon_name,
         return_giottoPolygon = FALSE
diff --git a/R/kriging.R b/R/kriging.R
index 2c0c8b204..95612a41f 100644
--- a/R/kriging.R
+++ b/R/kriging.R
@@ -182,8 +182,8 @@ setMethod(
         e_numeric <- ext[]
 
 
-        progressr::with_progress({
-            pb <- progressr::progressor(along = feats)
+        with_pbar({
+            pb <- pbar(along = feats)
 
             interp_img_list <- lapply_flex(
                 feats,
diff --git a/R/python_hmrf.R b/R/python_hmrf.R
index 6b0b597d0..a4ac65ed8 100644
--- a/R/python_hmrf.R
+++ b/R/python_hmrf.R
@@ -1498,7 +1498,7 @@ initHMRF_V2 <-
                 "scaled",
                 "normalized", "custom", expression_values
             )))
-            expr_values <- get_expression_values(
+            expr_values <- getExpression(
                 gobject = gobject,
                 spat_unit = spat_unit,
                 feat_type = feat_type,
@@ -1506,7 +1506,7 @@ initHMRF_V2 <-
             )
             if (zscore != "none") {
                 zscore <- match.arg(zscore, c("none", "colrow", "rowcol"))
-                expr_values <- get_expression_values(
+                expr_values <- getExpression(
                     gobject = gobject,
                     spat_unit = spat_unit,
                     feat_type = feat_type,
@@ -2034,7 +2034,7 @@ addHMRF_V2 <- function(gobject, HMRFoutput, name = "hmrf") {
             spat_unit = spat_unit,
             feat_type = feat_type,
             new_metadata = HMRFoutput[[i]]$class[match(
-                ordered_cell_IDs, 
+                ordered_cell_IDs,
                 rownames(HMRFoutput[[i]]$prob))],
             vector_name = paste(name, names(HMRFoutput)[i])
             # ,column_cell_ID = 'cell_ID',
@@ -2044,8 +2044,8 @@ addHMRF_V2 <- function(gobject, HMRFoutput, name = "hmrf") {
     return(gobject)
 }
 
-                             
-                             
+
+
 #' @title viewHMRFresults_V2
 #' @name viewHMRFresults_V2
 #' @description function to view HMRF results with multiple betas
@@ -2151,7 +2151,7 @@ viewHMRFresults_V2 <-
         }
 
         # combine plots with cowplot
-        combo_plot <- cowplot::plot_grid(
+        combo_plot <- plot_grid(
             plotlist = savelist,
             ncol = cow_n_col,
             rel_heights = cow_rel_h,
diff --git a/R/spatial_enrichment.R b/R/spatial_enrichment.R
index 1f95d9a53..5c92d38a0 100644
--- a/R/spatial_enrichment.R
+++ b/R/spatial_enrichment.R
@@ -1,38 +1,84 @@
 ## create spatial enrichment matrix ####
 
-#' @title makeSignMatrixPAGE
-#' @description Function to convert a list of signature genes
-#' (e.g. for cell types or processes) into
-#' a binary matrix format that can be used with the PAGE enrichment option.
-#' Each cell type or process should
-#' have a vector of cell-type or process specific genes. These vectors need to
-#' be combined into a list (sign_list).
-#' The names of the cell types or processes that are provided in the list need
-#' to be given (sign_names).
-#' @param sign_names vector with names for each provided gene signature
-#' @param sign_list list of genes (signature)
-#' @returns matrix
-#' @seealso \code{\link{PAGEEnrich}}
+#' @title PAGE feature enrichment
+#' @name enrichment_PAGE
+#' @description
+#' Expression feature-based enrichment scoring of labels. \cr
+#' A binary matrix of signature features (e.g. for cell types or processes) can
+#' either be directly provided or converted from a list using
+#' `makeSignMatrixPAGE()`. This matrix is then used with `runPAGEEnrich()` in
+#' order to calculate feature signature enrichment scores per spatial position
+#' using PAGE.
+#' @param sign_names `character` vector with names (labels) for each provided
+#' feat signature
+#' @param sign_list list of feats in signature
+#' @param gobject Giotto object
+#' @param spat_unit spatial unit
+#' @param feat_type feature type
+#' @param sign_matrix binary matrix of signature feats for each cell type /
+#' process. Alternatively a list of signature feats can be provided to
+#' `makeSignMatrixPAGE()`, which will create the matrix for you.
+#' @param expression_values expression values to use
+#' @param min_overlap_genes minimum number of overlapping feats in
+#' `sign_matrix` required to calculate enrichment
+#' @param reverse_log_scale reverse expression values from log scale
+#' @param logbase log base to use if reverse_log_scale = TRUE
+#' @param output_enrichment how to return enrichment output
+#' @param p_value logical. Default = `FALSE`. calculate p-values
+#' @param include_depletion calculate both enrichment and depletion
+#' @param n_times number of permutations to calculate for p_value
+#' @param max_block number of lines to process together (default = 20e6)
+#' @param name to give to spatial enrichment results, default = PAGE
+#' @param verbose be verbose
+#' @param return_gobject return giotto object
+#' @returns `matrix` (`makeSignMatrixPAGE()`) and
+#' `giotto` (`runPAGEEnrich(return_gobject = TRUE)`) or
+#' `data.table` (`runPAGEEnrich(return_gobject = FALSE)`)
+#' @details
+#' The enrichment Z score is calculated by using method (PAGE) from
+#' Kim SY et al., BMC bioinformatics, 2005 as \cr
+#' \eqn{Z = ((Sm – mu)*m^(1/2)) / delta}. \cr
+#' For each gene in each spot, mu is the fold change values versus the mean
+#' expression and delta is the standard deviation. Sm is the mean fold change
+#' value of a specific marker gene set and  m is the size of a given marker
+#' gene set.
+#' @md
 #' @examples
+#' g <- GiottoData::loadGiottoMini("visium")
+#'
 #' sign_list <- list(
 #'     cell_type1 = c(
 #'         "Bcl11b", "Lmo1", "F3", "Cnih3", "Ppp1r3c",
 #'         "Rims2", "Gfap", "Gjc3", "Chrna4", "Prkcd"
 #'     ),
 #'     cell_type2 = c(
-#'         "Prr18", "Grb14", "Tprn", "Clic1", "Olig2", "Hrh3", "Tmbim1",
-#'         "Carhsp1", "Tmem88b", "Ugt8a"
+#'         "Prr18", "Grb14", "Tprn", "Clic1", "Olig2", "Hrh3",
+#'          "Tmbim1", "Carhsp1", "Tmem88b", "Ugt8a"
 #'     ),
 #'     cell_type2 = c(
-#'         "Arpp19", "Lamp5", "Galnt6", "Hlf", "Hs3st2", "Tbr1", "Myl4",
-#'         "Cygb", "Ttc9b", "Ipcef1"
+#'         "Arpp19", "Lamp5", "Galnt6", "Hlf", "Hs3st2",
+#'         "Tbr1", "Myl4", "Cygb", "Ttc9b", "Ipcef1"
 #'     )
 #' )
 #'
-#' makeSignMatrixPAGE(
+#' sm <- makeSignMatrixPAGE(
 #'     sign_names = c("cell_type1", "cell_type2", "cell_type3"),
 #'     sign_list = sign_list
 #' )
+#'
+#' g <- runPAGEEnrich(gobject = g,
+#'     sign_matrix = sm,
+#'     min_overlap_genes = 2
+#' )
+#'
+#' spatPlot2D(g,
+#'     cell_color = "cell_type2",
+#'     spat_enr_names = "PAGE",
+#'     color_as_factor = FALSE
+#' )
+NULL
+
+#' @rdname enrichment_PAGE
 #' @export
 makeSignMatrixPAGE <- function(
         sign_names,
@@ -210,7 +256,7 @@ makeSignMatrixDWLS <- function(
         expression_values,
         unique(c("normalized", "scaled", "custom", expression_values))
     )
-    expr_values <- get_expression_values(
+    expr_values <- getExpression(
         gobject = gobject,
         spat_unit = spat_unit,
         feat_type = feat_type,
@@ -420,218 +466,6 @@ makeSignMatrixRank <- function(
 
 
 
-#' @title runPAGEEnrich_OLD
-#' @description Function to calculate gene signature enrichment scores per
-#' spatial position using PAGE.
-#' @param gobject Giotto object
-#' @param sign_matrix Matrix of signature genes for each cell type / process
-#' @param expression_values expression values to use
-#' @param reverse_log_scale reverse expression values from log scale
-#' @param logbase log base to use if reverse_log_scale = TRUE
-#' @param output_enrichment how to return enrichment output
-#' @param p_value calculate p-values (boolean, default = FALSE)
-#' @param n_times number of permutations to calculate for p_value
-#' @param name to give to spatial enrichment results, default = PAGE
-#' @param return_gobject return giotto object
-#' @returns data.table with enrichment results
-#' @details
-#' sign_matrix: a binary matrix with genes as row names and cell-types as
-#' column names.
-#' Alternatively a list of signature genes can be provided to
-#' makeSignMatrixPAGE, which will create the matrix for you. \cr
-#'
-#' The enrichment Z score is calculated by using method (PAGE) from
-#' Kim SY et al., BMC bioinformatics, 2005 as
-#' \eqn{Z = ((Sm – mu)*m^(1/2)) / delta}.
-#' For each gene in each spot, mu is the fold change values versus the mean
-#' expression and delta is the standard deviation. Sm is the mean fold change
-#' value of a specific marker gene set and  m is the size of a given marker
-#' gene set.
-#' @seealso \code{\link{makeSignMatrixPAGE}}
-#' @export
-runPAGEEnrich_OLD <- function(
-        gobject,
-        sign_matrix,
-        expression_values = c("normalized", "scaled", "custom"),
-        reverse_log_scale = TRUE,
-        logbase = 2,
-        output_enrichment = c("original", "zscore"),
-        p_value = FALSE,
-        n_times = 1000,
-        name = NULL,
-        return_gobject = TRUE) {
-    # expression values to be used
-    values <- match.arg(expression_values, c("normalized", "scaled", "custom"))
-    expr_values <- get_expression_values(gobject = gobject, values = values)
-
-    # check parameters
-    if (is.null(name)) name <- "PAGE"
-
-    # check available gene
-    available_ct <- c()
-    for (i in colnames(sign_matrix)) {
-        gene_i <- rownames(sign_matrix)[which(sign_matrix[, i] == 1)]
-        overlap_i <- intersect(gene_i, rownames(expr_values))
-        if (length(overlap_i) <= 5) {
-            output <- paste0(
-                "Warning, ", i, " only has ", length(overlap_i),
-                " overlapped genes. Will remove it."
-            )
-        } else {
-            available_ct <- c(available_ct, i)
-        }
-    }
-
-    if (length(available_ct) == 1) {
-        stop("Only one cell type available.")
-    }
-
-    # output enrichment
-    output_enrichment <- match.arg(
-        output_enrichment,
-        choices = c("original", "zscore")
-    )
-
-    # only continue with genes present in both datasets
-    interGene <- intersect(rownames(sign_matrix), rownames(expr_values))
-    filterSig <- sign_matrix[interGene, available_ct]
-    signames <- rownames(filterSig)[which(filterSig[, 1] == 1)]
-
-    # calculate mean gene expression
-    if (reverse_log_scale == TRUE) {
-        mean_gene_expr <- log(rowMeans(logbase^expr_values - 1, dims = 1) + 1)
-    } else {
-        mean_gene_expr <- rowMeans(expr_values)
-    }
-    geneFold <- expr_values - mean_gene_expr
-
-    # calculate sample/spot mean and sd
-    cellColMean <- apply(geneFold, 2, mean)
-    cellColSd <- apply(geneFold, 2, stats::sd)
-
-    # get enrichment scores
-    enrichment <- matrix(
-        data = NA, nrow = dim(filterSig)[2], ncol = length(cellColMean)
-    )
-    for (i in seq_len(dim(filterSig)[2])) {
-        signames <- rownames(filterSig)[which(filterSig[, i] == 1)]
-        sigColMean <- apply(geneFold[signames, ], 2, mean)
-        m <- length(signames)
-        vectorX <- NULL
-        for (j in (seq_along(cellColMean))) {
-            Sm <- sigColMean[j]
-            u <- cellColMean[j]
-            sigma <- cellColSd[j]
-            zscore <- (Sm - u) * m^(1 / 2) / sigma
-            vectorX <- append(vectorX, zscore)
-        }
-        enrichment[i, ] <- vectorX
-    }
-
-    rownames(enrichment) <- colnames(filterSig)
-    colnames(enrichment) <- names(cellColMean)
-    enrichment <- t(enrichment)
-
-    if (output_enrichment == "zscore") {
-        enrichment <- scale(enrichment)
-    }
-
-    enrichmentDT <- data.table::data.table(cell_ID = rownames(enrichment))
-    enrichmentDT <- cbind(enrichmentDT, data.table::as.data.table(enrichment))
-
-
-
-    ## calculate p-values if requested
-    if (p_value == TRUE) {
-        # check available gene
-        available_ct <- c()
-        for (i in colnames(sign_matrix)) {
-            gene_i <- rownames(sign_matrix)[which(sign_matrix[, i] == 1)]
-            overlap_i <- intersect(
-                gene_i, rownames(gobject@expression$rna$normalized)
-            )
-
-            if (length(overlap_i) <= 5) {
-                output <- paste0(
-                    "Warning, ", i, " only has ", length(overlap_i),
-                    " overlapped genes. It will be removed."
-                )
-            } else {
-                available_ct <- c(available_ct, i)
-            }
-        }
-
-        if (length(available_ct) == 1) {
-            stop("Only one cell type available.")
-        }
-
-        # only continue with genes present in both datasets
-        interGene <- intersect(
-            rownames(sign_matrix), rownames(gobject@expression$rna$normalized)
-        )
-        filter_sign_matrix <- sign_matrix[interGene, available_ct]
-
-        background_mean_sd <- .do_page_permutation(
-            gobject = gobject,
-            sig_gene = filter_sign_matrix,
-            ntimes = n_times
-        )
-
-        for (i in seq_len(dim(filter_sign_matrix)[2])) {
-            length_gene <- length(which(filter_sign_matrix[, i] == 1))
-            join_gene_with_length <- paste("gene_num_", length_gene, sep = "")
-            mean_i <- as.numeric(as.character(
-                background_mean_sd[join_gene_with_length, ][[1]]
-            ))
-            sd_i <- as.numeric(as.character(
-                background_mean_sd[join_gene_with_length, ][[2]]
-            ))
-            j <- i + 1
-            enrichmentDT[[j]] <- stats::pnorm(
-                enrichmentDT[[j]],
-                mean = mean_i, sd = sd_i,
-                lower.tail = FALSE, log.p = FALSE
-            )
-        }
-    }
-
-
-
-    ## return object or results ##
-    if (return_gobject == TRUE) {
-        spenr_names <- names(gobject@spatial_enrichment)
-
-        if (name %in% spenr_names) {
-            cat(name, " has already been used, will be overwritten")
-        }
-
-        ## update parameters used ##
-        parameters_list <- gobject@parameters
-        number_of_rounds <- length(parameters_list)
-        update_name <- paste0(number_of_rounds, "_spatial_enrichment")
-
-        # parameters to include
-        parameters_list[[update_name]] <- c(
-            "method used" = "PAGE",
-            "enrichment name" = name,
-            "expression values" = expression_values,
-            "reverse log scale" = reverse_log_scale,
-            "logbase" = logbase,
-            "p-values calculated" = p_value,
-            "output enrichment scores" = output_enrichment,
-            "p values calculated" = p_value,
-            "nr permutations" = n_times
-        )
-        gobject@parameters <- parameters_list
-
-        gobject@spatial_enrichment[[name]] <- enrichmentDT
-
-        return(gobject)
-    } else {
-        return(enrichmentDT)
-    }
-}
-
 
 
 #' @title PAGE data.table method
@@ -664,7 +498,7 @@ runPAGEEnrich_OLD <- function(
     ## identify available cell types
     all_genes <- rownames(expr_values)
     sign_matrix <- as.matrix(sign_matrix)
-    sign_matrix_DT <- data.table::as.data.table(reshape2::melt(sign_matrix))
+    sign_matrix_DT <- melt_matrix(sign_matrix)
     sign_matrix_DT <- sign_matrix_DT[Var1 %in% all_genes]
     detected_DT <- sign_matrix_DT[, sum(value), by = Var2]
 
@@ -685,7 +519,7 @@ runPAGEEnrich_OLD <- function(
         stop("Only one cell type available.")
     }
 
-    # create subset of sinature matrix
+    # create subset of signature matrix
     interGene <- intersect(rownames(sign_matrix), rownames(expr_values))
     filterSig <- sign_matrix[interGene, available_ct]
 
@@ -707,13 +541,13 @@ runPAGEEnrich_OLD <- function(
         colSd = cellColSd
     )
 
-    filterSig_DT <- data.table::as.data.table(reshape2::melt(filterSig))
+    filterSig_DT <- melt_matrix(filterSig)
     colnames(filterSig_DT) <- c("gene", "cell_type", "marker")
     sub_ct_DT <- filterSig_DT[marker == 1]
     sub_ct_DT[, nr_markers := .N, by = cell_type]
 
     ## reshape gene fold-expression
-    geneFold_DT <- data.table::as.data.table(reshape2::melt(geneFold))
+    geneFold_DT <- melt_matrix(geneFold)
     colnames(geneFold_DT) <- c("gene", "cell_ID", "fc")
 
     mergetest <- data.table::merge.data.table(
@@ -882,57 +716,7 @@ runPAGEEnrich_OLD <- function(
 
 
 
-#' @title runPAGEEnrich
-#' @description Function to calculate gene signature enrichment scores per
-#' spatial position using PAGE.
-#' @param gobject Giotto object
-#' @param spat_unit spatial unit
-#' @param feat_type feature type
-#' @param sign_matrix Matrix of signature genes for each cell type / process
-#' @param expression_values expression values to use
-#' @param min_overlap_genes minimum number of overlapping genes in sign_matrix
-#' required to calculate enrichment
-#' @param reverse_log_scale reverse expression values from log scale
-#' @param logbase log base to use if reverse_log_scale = TRUE
-#' @param output_enrichment how to return enrichment output
-#' @param p_value calculate p-values (boolean, default = FALSE)
-#' @param include_depletion calculate both enrichment and depletion
-#' @param n_times number of permutations to calculate for p_value
-#' @param max_block number of lines to process together (default = 20e6)
-#' @param name to give to spatial enrichment results, default = PAGE
-#' @param verbose be verbose
-#' @param return_gobject return giotto object
-#' @returns data.table with enrichment results
-#' @details
-#' sign_matrix: a binary matrix with genes as row names and cell-types as
-#' column names.
-#' Alternatively a list of signature genes can be provided to
-#' makeSignMatrixPAGE, which will create the matrix for you. \cr
-#'
-#' The enrichment Z score is calculated by using method (PAGE) from
-#' Kim SY et al., BMC bioinformatics, 2005 as
-#' \eqn{Z = ((Sm – mu)*m^(1/2)) / delta}.
-#' For each gene in each spot, mu is the fold change values versus the mean
-#' expression and delta is the standard deviation. Sm is the mean fold change
-#' value of a specific marker gene set and  m is the size of a given marker
-#' gene set.
-#' @seealso \code{\link{makeSignMatrixPAGE}}
-#' @examples
-#' g <- GiottoData::loadGiottoMini("visium")
-#' sign_gene <- c(
-#'     "Bcl11b", "Lmo1", "F3", "Cnih3", "Ppp1r3c", "Rims2", "Gfap",
-#'     "Gjc3", "Chrna4", "Prkcd", "Prr18", "Grb14", "Tprn", "Clic1", "Olig2",
-#'     "Hrh3", "Tmbim1", "Carhsp1", "Tmem88b", "Ugt8a", "Arpp19", "Lamp5",
-#'     "Galnt6", "Hlf", "Hs3st2", "Tbr1", "Myl4", "Cygb", "Ttc9b", "Ipcef1"
-#' )
-#'
-#' sign_matrix <- matrix(rnorm(length(sign_gene) * 3, mean = 10),
-#'     nrow = length(sign_gene)
-#' )
-#' rownames(sign_matrix) <- sign_gene
-#' colnames(sign_matrix) <- c("cell_type1", "cell_type2", "cell_type3")
-#'
-#' runPAGEEnrich(gobject = g, sign_matrix = sign_matrix)
+#' @rdname enrichment_PAGE
 #' @export
 runPAGEEnrich <- function(
         gobject,
@@ -1044,10 +828,7 @@ runPAGEEnrich <- function(
         gobject@parameters <- parameters_list
 
         ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ###
-        gobject <- set_spatial_enrichment(
-            gobject = gobject,
-            spatenrichment = enrObj
-        )
+        gobject <- setGiotto(gobject, enrObj)
         ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ###
 
 
@@ -1064,21 +845,6 @@ runPAGEEnrich <- function(
 
 
 
-#' @title PAGEEnrich
-#' @description Function to calculate gene signature enrichment scores per
-#' spatial position using PAGE.
-#' @inheritDotParams runPAGEEnrich
-#' @seealso \code{\link{runPAGEEnrich}}
-#' @returns gene enrichment
-#' @export
-PAGEEnrich <- function(...) {
-    .Deprecated(new = "runPAGEEnrich")
-
-    runPAGEEnrich(...)
-}
-
-
-
 
 
 #' @title Rank permutation
@@ -1363,10 +1129,7 @@ runRankEnrich <- function(
         gobject@parameters <- parameters_list
 
         ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ###
-        gobject <- set_spatial_enrichment(
-            gobject = gobject,
-            spatenrichment = enrObj
-        )
+        gobject <- setGiotto(gobject, enrObj)
         ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ###
 
         return(gobject)
@@ -1377,20 +1140,6 @@ runRankEnrich <- function(
 
 
 
-#' @title rankEnrich
-#' @description Function to calculate gene signature enrichment scores per
-#' spatial position using a rank based approach.
-#' @inheritDotParams runRankEnrich
-#' @seealso \code{\link{runRankEnrich}}
-#' @returns gene enrichment
-#' @export
-rankEnrich <- function(...) {
-    .Deprecated(new = "runRankEnrich")
-
-    runRankEnrich(...)
-}
-
-
 
 #' @title runHyperGeometricEnrich
 #' @description Function to calculate gene signature enrichment scores per
@@ -1596,10 +1345,7 @@ runHyperGeometricEnrich <- function(
         gobject@parameters <- parameters_list
 
         ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ###
-        gobject <- set_spatial_enrichment(
-            gobject = gobject,
-            spatenrichment = enrObj
-        )
+        gobject <- setGiotto(gobject, enrObj)
         ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ###
 
         return(gobject)
@@ -1609,19 +1355,6 @@ runHyperGeometricEnrich <- function(
 }
 
 
-#' @title hyperGeometricEnrich
-#' @description Function to calculate gene signature enrichment scores per
-#' spatial position using a hypergeometric test.
-#' @inheritDotParams runHyperGeometricEnrich
-#' @seealso \code{\link{runHyperGeometricEnrich}}
-#' @returns enrichment scores
-#' @export
-hyperGeometricEnrich <- function(...) {
-    .Deprecated(new = "runHyperGeometricEnrich")
-
-    runHyperGeometricEnrich(...)
-}
-
 
 
 
@@ -1755,18 +1488,6 @@ runSpatialEnrich <- function(
 }
 
 
-#' @title createSpatialEnrich
-#' @description Function to calculate gene signature enrichment scores per
-#' spatial position using an enrichment test.
-#' @inheritDotParams runSpatialEnrich
-#' @returns gene signature enrichment scores
-#' @seealso \code{\link{runSpatialEnrich}}
-#' @export
-createSpatialEnrich <- function(...) {
-    .Deprecated(new = "runSpatialEnrich")
-
-    runSpatialEnrich(...)
-}
 
 
 
@@ -2174,10 +1895,7 @@ spatialAutoCorLocal <- function(
             )
         }
         ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ###
-        gobject <- set_spatial_enrichment(
-            gobject = gobject,
-            spatenrichment = enr
-        )
+        gobject <- setGiotto(gobject, enr)
         ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ###
 
         return(gobject)
@@ -2218,9 +1936,9 @@ spatialAutoCorLocal <- function(
         step_size <- step_size <- ceiling(nfeats / 10L)
     }
 
-    progressr::with_progress({
+    with_pbar({
         if (step_size > 1) {
-            pb <- progressr::progressor(
+            pb <- pbar(
                 steps = nfeats / step_size
             )
         }
@@ -2301,9 +2019,9 @@ spatialAutoCorLocal <- function(
         step_size <- step_size <- ceiling(nfeats / 10L)
     }
 
-    progressr::with_progress({
+    with_pbar({
         if (step_size > 1) {
-            pb <- progressr::progressor(
+            pb <- pbar(
                 steps = nfeats / step_size
             )
         }
@@ -3198,10 +2916,7 @@ runDWLSDeconv <- function(
         gobject@parameters <- parameters_list
 
         ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ###
-        gobject <- set_spatial_enrichment(
-            gobject = gobject,
-            spatenrichment = enrObj
-        )
+        gobject <- setGiotto(gobject, enrObj)
         ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ###
 
         return(gobject)
diff --git a/R/spatial_genes.R b/R/spatial_genes.R
index f77cc1527..473ecaadc 100644
--- a/R/spatial_genes.R
+++ b/R/spatial_genes.R
@@ -4104,19 +4104,6 @@ heatmSpatialCorFeats <- function(
 
 
 
-#' @title heatmSpatialCorGenes
-#' @name heatmSpatialCorGenes
-#' @description Create heatmap of spatially correlated genes
-#' @inheritDotParams heatmSpatialCorFeats
-#' @returns heatmap
-#' @seealso \code{\link{heatmSpatialCorFeats}}
-#' @export
-heatmSpatialCorGenes <- function(...) {
-    .Deprecated(new = "heatmSpatialCorFeats")
-
-    heatmSpatialCorFeats(...)
-}
-
 
 
 
diff --git a/R/spatial_interaction.R b/R/spatial_interaction.R
index 27e82edea..c009c84a2 100644
--- a/R/spatial_interaction.R
+++ b/R/spatial_interaction.R
@@ -1996,8 +1996,10 @@ print.combIcfObject <- function(x, ...) {
         lig_cell_type <- rec_cell_type <- NULL
 
     all_ligand_cols <- colnames(ligand_match)
-    lig_test <- data.table::as.data.table(
-        reshape2::melt(ligand_match, measure.vars = all_ligand_cols))
+    lig_test <- data.table::melt(
+        data.table::as.data.table(ligand_match),
+        measure.vars = all_ligand_cols
+    )
     lig_test[, ligand := rep(rownames(ligand_match), ncol(ligand_match))]
     lig_test[, ligand := strsplit(ligand, "\\.")[[1]][1],
             by = seq_len(nrow(lig_test))]
@@ -2006,8 +2008,10 @@ print.combIcfObject <- function(x, ...) {
     setnames(lig_test, "variable", "lig_cell_type")
 
     all_receptor_cols <- colnames(receptor_match)
-    rec_test <- data.table::as.data.table(reshape2::melt(
-        receptor_match, measure.vars = all_receptor_cols))
+    rec_test <- data.table::melt(
+        data.table::as.data.table(receptor_match),
+        measure.vars = all_receptor_cols
+    )
     rec_test[, receptor := rep(rownames(receptor_match), ncol(receptor_match))]
     rec_test[, receptor := strsplit(
         receptor, "\\.")[[1]][1], by = seq_len(nrow(rec_test))]
@@ -2045,7 +2049,7 @@ print.combIcfObject <- function(x, ...) {
 #' (random variance and z-score)
 #' @param adjust_method which method to adjust p-values
 #' @param adjust_target adjust multiple hypotheses at the cell or feature level
-#' @param set_seed set seed for random simulations (default = TRUE)
+#' @param set_seed `logical`. set seed for random simulations (default = TRUE)
 #' @param seed_number seed number
 #' @param verbose verbose
 #' @returns Cell-Cell communication scores for feature pairs based on
@@ -2145,8 +2149,8 @@ exprCellCellcom <- function(gobject,
     # not yet available
 
 
-    progressr::with_progress({
-        pb <- progressr::progressor(steps = random_iter)
+    with_pbar({
+        pb <- pbar(steps = random_iter)
 
         for (sim in seq_len(random_iter)) {
 
diff --git a/R/spatial_interaction_spot.R b/R/spatial_interaction_spot.R
index ab8ed82d4..267c426c8 100644
--- a/R/spatial_interaction_spot.R
+++ b/R/spatial_interaction_spot.R
@@ -54,7 +54,7 @@ NULL
             diag(diff_ct) <- NA
             diff_ct[lower.tri(diff_ct)] <- NA
             # transfer format to data.table
-            diff_ct <- data.table::as.data.table(reshape2::melt(diff_ct))
+            diff_ct <- melt_matrix(diff_ct)
             diff_ct <- diff_ct[value != "NA"]
             diff_ct[, c("Var1", "Var2") := lapply(
                 .SD, as.character
@@ -99,10 +99,9 @@ NULL
     pairs_for_mat <- pairs_for_mat[, .N, by = c("from", "to")]
 
     # make square matrix of interaction between spots
-    pairs_mat <- reshape2::acast(
-        pairs_for_mat, from ~ to,
-        value.var = "N", fill = 0
-    )
+    pairs_mat <- as.matrix(data.table::dcast(
+        pairs_for_mat, from ~ to, value.var = "N", fill = 0
+    ), rownames = "from")
     pairs_mat <- pairs_mat[cell_IDs, cell_IDs]
 
     # calculate cell-type/cell-type interactions
@@ -225,7 +224,7 @@ NULL
 #' )
 #' sign_gene <- x$feats
 #'
-#' sign_matrix <- matrix(rnorm(length(sign_gene) * 8, mean = 10),
+#' sign_matrix <- matrix(rnorm(length(sign_gene) * 7, mean = 10),
 #'     nrow = length(sign_gene)
 #' )
 #' rownames(sign_matrix) <- sign_gene
@@ -715,8 +714,8 @@ cellProximityEnrichmentEachSpot <- function(
                 dimnames = list(names(dwls_target_cell), names(spot_proximity))
             )
         }
-        spot_proximity <- reshape2::melt(spot_proximity)
-        spot_proximity <- data.table::data.table(spot_proximity)
+
+        spot_proximity <- melt_matrix(spot_proximity)
         spot_proximity[, c("Var1", "Var2") := lapply(
             .SD, as.character
         ), .SDcols = c("Var1", "Var2")]
@@ -1146,7 +1145,9 @@ NULL
 #' @param spat_unit spatial unit (e.g. 'cell')
 #' @param feat_type feature type (e.g. 'rna')
 #' @param expression_values expression values to use
-#' @param ave_celltype_exp average feature expression in each cell type
+#' @param ave_celltype_exp `matrix` or `data.frame`. Average feature expression
+#' in each cell type. Colnames should be the cell type and rownames are feat
+#' names.
 #' @param selected_features subset of selected features (optional)
 #' @param spatial_network_name name of spatial network to use
 #' @param deconv_name name of deconvolution/spatial enrichment values to use
@@ -1197,7 +1198,7 @@ NULL
 #' )
 #' sign_gene <- x$feats
 #'
-#' sign_matrix <- matrix(rnorm(length(sign_gene) * 8, mean = 10),
+#' sign_matrix <- matrix(rnorm(length(sign_gene) * 7, mean = 10),
 #'     nrow = length(sign_gene)
 #' )
 #' rownames(sign_matrix) <- sign_gene
@@ -1205,17 +1206,19 @@ NULL
 #'
 #' g <- runDWLSDeconv(gobject = g, sign_matrix = sign_matrix)
 #' ave_celltype_exp <- calculateMetaTable(g, metadata_cols = "leiden_clus")
-#' ave_celltype_exp <- reshape2::dcast(ave_celltype_exp, variable~leiden_clus)
-#' rownames(ave_celltype_exp) <- ave_celltype_exp$variable
-#' ave_celltype_exp <- ave_celltype_exp[,-1]
+#' ave_celltype_exp <- data.table::dcast(
+#'     ave_celltype_exp, variable~leiden_clus
+#' )
+#' ave_celltype_exp <- as.matrix(ave_celltype_exp, rownames = "variable")
 #' colnames(ave_celltype_exp) <- colnames(sign_matrix)
-#' 
-#' findICFSpot(g,
+#'
+#' res <- findICFSpot(g,
 #'     spat_unit = "cell",
 #'     feat_type = "rna",
 #'     ave_celltype_exp = ave_celltype_exp,
 #'     spatial_network_name = "spatial_network"
 #' )
+#' force(res)
 #' @seealso [findInteractionChangedFeats()]
 #' @md
 #' @export
@@ -2252,6 +2255,8 @@ plotCellProximityFeatSpot <- function(
 #' @param expression_values (e.g. 'normalized', 'scaled', 'custom')
 #' @param spatial_network_name spatial network to use for identifying
 #' interacting cells
+#' @param spat_enr_name name of spatial enrichment containing DWLS results.
+#' Default = `"DWLS"`
 #' @param cluster_column cluster column with cell type information
 #' @param random_iter number of iterations
 #' @param feature_set_1 first specific feature set from feature pairs
@@ -2303,6 +2308,7 @@ plotCellProximityFeatSpot <- function(
 #'  * p.adj: adjusted p-value
 #'  * PI: significanc score: log2fc \* -log10(p.adj)
 #' }
+#' @md
 #' @export
 spatCellCellcomSpots <- function(
         gobject,
@@ -2310,6 +2316,7 @@ spatCellCellcomSpots <- function(
         feat_type = NULL,
         ave_celltype_exp,
         spatial_network_name = "Delaunay_network",
+        spat_enr_name = "DWLS",
         cluster_column = "cell_ID",
         random_iter = 1000,
         feature_set_1,
@@ -2378,10 +2385,11 @@ spatCellCellcomSpots <- function(
     proximityMat <- proximityMat[, intersect_cell_IDs]
 
     # exact spatial_enrichment matrix
-    dwls_values <- get_spatial_enrichment(
+    dwls_values <- getSpatialEnrichment(
         gobject = gobject,
         spat_unit = spat_unit,
         feat_type = feat_type,
+        name = spat_enr_name,
         output = "data.table"
     )
     data.table::setDF(dwls_values)
diff --git a/R/spatial_interaction_visuals.R b/R/spatial_interaction_visuals.R
index c74da4fc3..8f192b778 100644
--- a/R/spatial_interaction_visuals.R
+++ b/R/spatial_interaction_visuals.R
@@ -68,7 +68,7 @@ cellProximityBarplot <- function(
     bpl <- bpl + ggplot2::labs(y = "# of interactions")
     bpl
 
-    combo_plot <- cowplot::plot_grid(
+    combo_plot <- plot_grid(
         pl, bpl,
         ncol = 2, rel_heights = c(1),
         rel_widths = c(3, 1.5), align = "h"
diff --git a/R/spdep.R b/R/spdep.R
index 6abbe9cad..84b1b9ed9 100644
--- a/R/spdep.R
+++ b/R/spdep.R
@@ -77,9 +77,9 @@ spdepAutoCorr <- function(
     step_size <- ceiling(nfeats / 10L)
 
     result_list <- list()
-    progressr::with_progress({
+    with_pbar({
         if (step_size > 1) {
-            pb <- progressr::progressor(
+            pb <- pbar(
                 steps = nfeats / step_size
             )
         }
diff --git a/R/wnn.R b/R/wnn.R
index 588d11888..25b2cb829 100644
--- a/R/wnn.R
+++ b/R/wnn.R
@@ -458,7 +458,7 @@ runWNN <- function(
         matrix_result_name <- "theta_weighted_matrix"
     }
 
-    gobject <- set_multiomics(
+    gobject <- setMultiomics(
         gobject = gobject,
         result = theta_weighted,
         spat_unit = spat_unit,
@@ -476,7 +476,7 @@ runWNN <- function(
         w_name_modality_1 <- paste0("w_", modality_1)
     }
 
-    gobject <- set_multiomics(
+    gobject <- setMultiomics(
         gobject = gobject,
         result = w_modality1,
         spat_unit = spat_unit,
@@ -491,7 +491,7 @@ runWNN <- function(
         w_name_modality_2 <- paste0("w_", modality_2)
     }
 
-    gobject <- set_multiomics(
+    gobject <- setMultiomics(
         gobject = gobject,
         result = w_modality2,
         spat_unit = spat_unit,
@@ -539,7 +539,7 @@ runIntegratedUMAP <- function(
         integrated_feat_type <- paste0(modality1, "_", modality2)
     }
 
-    theta_weighted <- get_multiomics(gobject,
+    theta_weighted <- getMultiomics(gobject,
         spat_unit = spat_unit,
         feat_type = integrated_feat_type,
         integration_method = integration_method,
@@ -594,17 +594,10 @@ runIntegratedUMAP <- function(
             feat_type = modality1
         )
 
-        gobject <- set_NearestNetwork(
-            gobject = gobject,
-            nn_network = nnNetObj,
-            spat_unit = spat_unit,
-            feat_type = modality1,
-            nn_network_to_use = "kNN",
-            network_name = "integrated_kNN"
-        )
+        gobject <- setGiotto(gobject, nnNetObj)
 
         ## store nn_network id
-        gobject <- set_multiomics(
+        gobject <- setMultiomics(
             gobject = gobject,
             result = nn_network$id,
             spat_unit = spat_unit,
@@ -615,7 +608,7 @@ runIntegratedUMAP <- function(
         )
 
         ## store nn_network dist
-        gobject <- set_multiomics(
+        gobject <- setMultiomics(
             gobject = gobject,
             result = nn_network$dist,
             spat_unit = spat_unit,
diff --git a/man/PAGEEnrich.Rd b/man/PAGEEnrich.Rd
deleted file mode 100644
index 50596ada6..000000000
--- a/man/PAGEEnrich.Rd
+++ /dev/null
@@ -1,41 +0,0 @@
-% Generated by roxygen2: do not edit by hand
-% Please edit documentation in R/spatial_enrichment.R
-\name{PAGEEnrich}
-\alias{PAGEEnrich}
-\title{PAGEEnrich}
-\usage{
-PAGEEnrich(...)
-}
-\arguments{
-\item{...}{
-  Arguments passed on to \code{\link[=runPAGEEnrich]{runPAGEEnrich}}
-  \describe{
-    \item{\code{gobject}}{Giotto object}
-    \item{\code{spat_unit}}{spatial unit}
-    \item{\code{feat_type}}{feature type}
-    \item{\code{sign_matrix}}{Matrix of signature genes for each cell type / process}
-    \item{\code{expression_values}}{expression values to use}
-    \item{\code{min_overlap_genes}}{minimum number of overlapping genes in sign_matrix
-required to calculate enrichment}
-    \item{\code{reverse_log_scale}}{reverse expression values from log scale}
-    \item{\code{logbase}}{log base to use if reverse_log_scale = TRUE}
-    \item{\code{output_enrichment}}{how to return enrichment output}
-    \item{\code{p_value}}{calculate p-values (boolean, default = FALSE)}
-    \item{\code{include_depletion}}{calculate both enrichment and depletion}
-    \item{\code{n_times}}{number of permutations to calculate for p_value}
-    \item{\code{max_block}}{number of lines to process together (default = 20e6)}
-    \item{\code{name}}{to give to spatial enrichment results, default = PAGE}
-    \item{\code{verbose}}{be verbose}
-    \item{\code{return_gobject}}{return giotto object}
-  }}
-}
-\value{
-gene enrichment
-}
-\description{
-Function to calculate gene signature enrichment scores per
-spatial position using PAGE.
-}
-\seealso{
-\code{\link{runPAGEEnrich}}
-}
diff --git a/man/cellProximityEnrichmentSpots.Rd b/man/cellProximityEnrichmentSpots.Rd
index e618dbbc3..432921401 100644
--- a/man/cellProximityEnrichmentSpots.Rd
+++ b/man/cellProximityEnrichmentSpots.Rd
@@ -70,7 +70,7 @@ x <- findMarkers_one_vs_all(g,
 )
 sign_gene <- x$feats
 
-sign_matrix <- matrix(rnorm(length(sign_gene) * 8, mean = 10),
+sign_matrix <- matrix(rnorm(length(sign_gene) * 7, mean = 10),
     nrow = length(sign_gene)
 )
 rownames(sign_matrix) <- sign_gene
diff --git a/man/createSpatialEnrich.Rd b/man/createSpatialEnrich.Rd
deleted file mode 100644
index c62e820da..000000000
--- a/man/createSpatialEnrich.Rd
+++ /dev/null
@@ -1,46 +0,0 @@
-% Generated by roxygen2: do not edit by hand
-% Please edit documentation in R/spatial_enrichment.R
-\name{createSpatialEnrich}
-\alias{createSpatialEnrich}
-\title{createSpatialEnrich}
-\usage{
-createSpatialEnrich(...)
-}
-\arguments{
-\item{...}{
-  Arguments passed on to \code{\link[=runSpatialEnrich]{runSpatialEnrich}}
-  \describe{
-    \item{\code{gobject}}{Giotto object}
-    \item{\code{spat_unit}}{spatial unit}
-    \item{\code{feat_type}}{feature type}
-    \item{\code{enrich_method}}{method for gene signature enrichment calculation}
-    \item{\code{sign_matrix}}{Matrix of signature genes for each cell type / process}
-    \item{\code{expression_values}}{expression values to use}
-    \item{\code{reverse_log_scale}}{reverse expression values from log scale}
-    \item{\code{min_overlap_genes}}{minimum number of overlapping genes in sign_matrix
-required to calculate enrichment (PAGE)}
-    \item{\code{logbase}}{log base to use if reverse_log_scale = TRUE}
-    \item{\code{p_value}}{calculate p-value (default = FALSE)}
-    \item{\code{n_times}}{(page/rank) number of permutation iterations to calculate
-p-value}
-    \item{\code{rbp_p}}{(rank) fractional binarization threshold (default = 0.99)}
-    \item{\code{num_agg}}{(rank) number of top genes to aggregate (default = 100)}
-    \item{\code{max_block}}{number of lines to process together (default = 20e6)}
-    \item{\code{top_percentage}}{(hyper) percentage of cells that will be considered
-to have gene expression with matrix binarization}
-    \item{\code{output_enrichment}}{how to return enrichment output}
-    \item{\code{name}}{to give to spatial enrichment results, default = PAGE}
-    \item{\code{verbose}}{be verbose}
-    \item{\code{return_gobject}}{return giotto object}
-  }}
-}
-\value{
-gene signature enrichment scores
-}
-\description{
-Function to calculate gene signature enrichment scores per
-spatial position using an enrichment test.
-}
-\seealso{
-\code{\link{runSpatialEnrich}}
-}
diff --git a/man/doKmeans.Rd b/man/doKmeans.Rd
index 3bfc8552d..ca9fde554 100644
--- a/man/doKmeans.Rd
+++ b/man/doKmeans.Rd
@@ -71,8 +71,8 @@ if return_gobject = TRUE: giotto object with new clusters appended to cell metad
 cluster cells using kmeans algorithm
 }
 \details{
-The default settings will use dimension reduction results as input. 
-Set dim_reduction_to_use = NULL if you want to directly use expression values as input. 
+The default settings will use dimension reduction results as input.
+Set dim_reduction_to_use = NULL if you want to directly use expression values as input.
 By providing a feature vector to feats_to_use you can subset the expression matrix.
 }
 \examples{
diff --git a/man/enrichment_PAGE.Rd b/man/enrichment_PAGE.Rd
new file mode 100644
index 000000000..226fcb31a
--- /dev/null
+++ b/man/enrichment_PAGE.Rd
@@ -0,0 +1,126 @@
+% Generated by roxygen2: do not edit by hand
+% Please edit documentation in R/spatial_enrichment.R
+\name{enrichment_PAGE}
+\alias{enrichment_PAGE}
+\alias{makeSignMatrixPAGE}
+\alias{runPAGEEnrich}
+\title{PAGE feature enrichment}
+\usage{
+makeSignMatrixPAGE(sign_names, sign_list)
+
+runPAGEEnrich(
+  gobject,
+  spat_unit = NULL,
+  feat_type = NULL,
+  sign_matrix,
+  expression_values = c("normalized", "scaled", "custom"),
+  min_overlap_genes = 5,
+  reverse_log_scale = TRUE,
+  logbase = 2,
+  output_enrichment = c("original", "zscore"),
+  p_value = FALSE,
+  include_depletion = FALSE,
+  n_times = 1000,
+  max_block = 2e+07,
+  name = NULL,
+  verbose = TRUE,
+  return_gobject = TRUE
+)
+}
+\arguments{
+\item{sign_names}{\code{character} vector with names (labels) for each provided
+feat signature}
+
+\item{sign_list}{list of feats in signature}
+
+\item{gobject}{Giotto object}
+
+\item{spat_unit}{spatial unit}
+
+\item{feat_type}{feature type}
+
+\item{sign_matrix}{binary matrix of signature feats for each cell type /
+process. Alternatively a list of signature feats can be provided to
+\code{makeSignMatrixPAGE()}, which will create the matrix for you.}
+
+\item{expression_values}{expression values to use}
+
+\item{min_overlap_genes}{minimum number of overlapping feats in
+\code{sign_matrix} required to calculate enrichment}
+
+\item{reverse_log_scale}{reverse expression values from log scale}
+
+\item{logbase}{log base to use if reverse_log_scale = TRUE}
+
+\item{output_enrichment}{how to return enrichment output}
+
+\item{p_value}{logical. Default = \code{FALSE}. calculate p-values}
+
+\item{include_depletion}{calculate both enrichment and depletion}
+
+\item{n_times}{number of permutations to calculate for p_value}
+
+\item{max_block}{number of lines to process together (default = 20e6)}
+
+\item{name}{to give to spatial enrichment results, default = PAGE}
+
+\item{verbose}{be verbose}
+
+\item{return_gobject}{return giotto object}
+}
+\value{
+\code{matrix} (\code{makeSignMatrixPAGE()}) and
+\code{giotto} (\code{runPAGEEnrich(return_gobject = TRUE)}) or
+\code{data.table} (\code{runPAGEEnrich(return_gobject = FALSE)})
+}
+\description{
+Expression feature-based enrichment scoring of labels. \cr
+A binary matrix of signature features (e.g. for cell types or processes) can
+either be directly provided or converted from a list using
+\code{makeSignMatrixPAGE()}. This matrix is then used with \code{runPAGEEnrich()} in
+order to calculate feature signature enrichment scores per spatial position
+using PAGE.
+}
+\details{
+The enrichment Z score is calculated by using method (PAGE) from
+Kim SY et al., BMC bioinformatics, 2005 as \cr
+\eqn{Z = ((Sm – mu)*m^(1/2)) / delta}. \cr
+For each gene in each spot, mu is the fold change values versus the mean
+expression and delta is the standard deviation. Sm is the mean fold change
+value of a specific marker gene set and  m is the size of a given marker
+gene set.
+}
+\examples{
+g <- GiottoData::loadGiottoMini("visium")
+
+sign_list <- list(
+    cell_type1 = c(
+        "Bcl11b", "Lmo1", "F3", "Cnih3", "Ppp1r3c",
+        "Rims2", "Gfap", "Gjc3", "Chrna4", "Prkcd"
+    ),
+    cell_type2 = c(
+        "Prr18", "Grb14", "Tprn", "Clic1", "Olig2", "Hrh3",
+         "Tmbim1", "Carhsp1", "Tmem88b", "Ugt8a"
+    ),
+    cell_type2 = c(
+        "Arpp19", "Lamp5", "Galnt6", "Hlf", "Hs3st2",
+        "Tbr1", "Myl4", "Cygb", "Ttc9b", "Ipcef1"
+    )
+)
+
+sm <- makeSignMatrixPAGE(
+    sign_names = c("cell_type1", "cell_type2", "cell_type3"),
+    sign_list = sign_list
+)
+
+g <- runPAGEEnrich(gobject = g,
+    sign_matrix = sm,
+    min_overlap_genes = 2
+)
+
+spatPlot2D(g,
+    cell_color = "cell_type2",
+    spat_enr_names = "PAGE",
+    color_as_factor = FALSE
+)
+}
diff --git a/man/exprCellCellcom.Rd b/man/exprCellCellcom.Rd
index e3d89cef8..ca167a8c1 100644
--- a/man/exprCellCellcom.Rd
+++ b/man/exprCellCellcom.Rd
@@ -46,7 +46,7 @@ exprCellCellcom(
 
 \item{adjust_target}{adjust multiple hypotheses at the cell or feature level}
 
-\item{set_seed}{set seed for random simulations (default = TRUE)}
+\item{set_seed}{`logical`. set seed for random simulations (default = TRUE)}
 
 \item{seed_number}{seed number}
 
diff --git a/man/findICFSpot.Rd b/man/findICFSpot.Rd
index 500f53e2d..91eeedf11 100644
--- a/man/findICFSpot.Rd
+++ b/man/findICFSpot.Rd
@@ -36,7 +36,9 @@ findICFSpot(
 
 \item{expression_values}{expression values to use}
 
-\item{ave_celltype_exp}{average feature expression in each cell type}
+\item{ave_celltype_exp}{\code{matrix} or \code{data.frame}. Average feature expression
+in each cell type. Colnames should be the cell type and rownames are feat
+names.}
 
 \item{selected_features}{subset of selected features (optional)}
 
@@ -113,7 +115,7 @@ x <- findMarkers_one_vs_all(g,
 )
 sign_gene <- x$feats
 
-sign_matrix <- matrix(rnorm(length(sign_gene) * 8, mean = 10),
+sign_matrix <- matrix(rnorm(length(sign_gene) * 7, mean = 10),
     nrow = length(sign_gene)
 )
 rownames(sign_matrix) <- sign_gene
@@ -121,17 +123,19 @@ colnames(sign_matrix) <- paste0("celltype_", unique(x$cluster))
 
 g <- runDWLSDeconv(gobject = g, sign_matrix = sign_matrix)
 ave_celltype_exp <- calculateMetaTable(g, metadata_cols = "leiden_clus")
-ave_celltype_exp <- reshape2::dcast(ave_celltype_exp, variable~leiden_clus)
-rownames(ave_celltype_exp) <- ave_celltype_exp$variable
-ave_celltype_exp <- ave_celltype_exp[,-1]
+ave_celltype_exp <- data.table::dcast(
+    ave_celltype_exp, variable~leiden_clus
+)
+ave_celltype_exp <- as.matrix(ave_celltype_exp, rownames = "variable")
 colnames(ave_celltype_exp) <- colnames(sign_matrix)
 
-findICFSpot(g,
+res <- findICFSpot(g,
     spat_unit = "cell",
     feat_type = "rna",
     ave_celltype_exp = ave_celltype_exp,
     spatial_network_name = "spatial_network"
 )
+force(res)
 }
 \seealso{
 \code{\link[=findInteractionChangedFeats]{findInteractionChangedFeats()}}
diff --git a/man/heatmSpatialCorGenes.Rd b/man/heatmSpatialCorGenes.Rd
deleted file mode 100644
index bbf314b0a..000000000
--- a/man/heatmSpatialCorGenes.Rd
+++ /dev/null
@@ -1,39 +0,0 @@
-% Generated by roxygen2: do not edit by hand
-% Please edit documentation in R/spatial_genes.R
-\name{heatmSpatialCorGenes}
-\alias{heatmSpatialCorGenes}
-\title{heatmSpatialCorGenes}
-\usage{
-heatmSpatialCorGenes(...)
-}
-\arguments{
-\item{...}{
-  Arguments passed on to \code{\link[=heatmSpatialCorFeats]{heatmSpatialCorFeats}}
-  \describe{
-    \item{\code{gobject}}{giotto object}
-    \item{\code{spatCorObject}}{spatial correlation object}
-    \item{\code{use_clus_name}}{name of clusters to visualize
-(from clusterSpatialCorFeats())}
-    \item{\code{show_cluster_annot}}{show cluster annotation on top of heatmap}
-    \item{\code{show_row_dend}}{show row dendrogram}
-    \item{\code{show_column_dend}}{show column dendrogram}
-    \item{\code{show_row_names}}{show row names}
-    \item{\code{show_column_names}}{show column names}
-    \item{\code{show_plot}}{show plot}
-    \item{\code{return_plot}}{return ggplot object}
-    \item{\code{save_plot}}{directly save the plot [boolean]}
-    \item{\code{save_param}}{list of saving parameters, see
-\code{\link{showSaveParameters}}}
-    \item{\code{default_save_name}}{default save name for saving, don't change,
-change save_name in save_param}
-  }}
-}
-\value{
-heatmap
-}
-\description{
-Create heatmap of spatially correlated genes
-}
-\seealso{
-\code{\link{heatmSpatialCorFeats}}
-}
diff --git a/man/hyperGeometricEnrich.Rd b/man/hyperGeometricEnrich.Rd
deleted file mode 100644
index c4dd29587..000000000
--- a/man/hyperGeometricEnrich.Rd
+++ /dev/null
@@ -1,37 +0,0 @@
-% Generated by roxygen2: do not edit by hand
-% Please edit documentation in R/spatial_enrichment.R
-\name{hyperGeometricEnrich}
-\alias{hyperGeometricEnrich}
-\title{hyperGeometricEnrich}
-\usage{
-hyperGeometricEnrich(...)
-}
-\arguments{
-\item{...}{
-  Arguments passed on to \code{\link[=runHyperGeometricEnrich]{runHyperGeometricEnrich}}
-  \describe{
-    \item{\code{gobject}}{Giotto object}
-    \item{\code{spat_unit}}{spatial unit}
-    \item{\code{feat_type}}{feature type}
-    \item{\code{sign_matrix}}{Matrix of signature genes for each cell type / process}
-    \item{\code{expression_values}}{expression values to use}
-    \item{\code{reverse_log_scale}}{reverse expression values from log scale}
-    \item{\code{logbase}}{log base to use if reverse_log_scale = TRUE}
-    \item{\code{top_percentage}}{percentage of cells that will be considered to have
-gene expression with matrix binarization}
-    \item{\code{output_enrichment}}{how to return enrichment output}
-    \item{\code{p_value}}{calculate p-values (boolean, default = FALSE)}
-    \item{\code{name}}{to give to spatial enrichment results, default = hypergeometric}
-    \item{\code{return_gobject}}{return giotto object}
-  }}
-}
-\value{
-enrichment scores
-}
-\description{
-Function to calculate gene signature enrichment scores per
-spatial position using a hypergeometric test.
-}
-\seealso{
-\code{\link{runHyperGeometricEnrich}}
-}
diff --git a/man/makeSignMatrixPAGE.Rd b/man/makeSignMatrixPAGE.Rd
deleted file mode 100644
index 805dde245..000000000
--- a/man/makeSignMatrixPAGE.Rd
+++ /dev/null
@@ -1,50 +0,0 @@
-% Generated by roxygen2: do not edit by hand
-% Please edit documentation in R/spatial_enrichment.R
-\name{makeSignMatrixPAGE}
-\alias{makeSignMatrixPAGE}
-\title{makeSignMatrixPAGE}
-\usage{
-makeSignMatrixPAGE(sign_names, sign_list)
-}
-\arguments{
-\item{sign_names}{vector with names for each provided gene signature}
-
-\item{sign_list}{list of genes (signature)}
-}
-\value{
-matrix
-}
-\description{
-Function to convert a list of signature genes
-(e.g. for cell types or processes) into
-a binary matrix format that can be used with the PAGE enrichment option.
-Each cell type or process should
-have a vector of cell-type or process specific genes. These vectors need to
-be combined into a list (sign_list).
-The names of the cell types or processes that are provided in the list need
-to be given (sign_names).
-}
-\examples{
-sign_list <- list(
-    cell_type1 = c(
-        "Bcl11b", "Lmo1", "F3", "Cnih3", "Ppp1r3c",
-        "Rims2", "Gfap", "Gjc3", "Chrna4", "Prkcd"
-    ),
-    cell_type2 = c(
-        "Prr18", "Grb14", "Tprn", "Clic1", "Olig2", "Hrh3", "Tmbim1",
-        "Carhsp1", "Tmem88b", "Ugt8a"
-    ),
-    cell_type2 = c(
-        "Arpp19", "Lamp5", "Galnt6", "Hlf", "Hs3st2", "Tbr1", "Myl4",
-        "Cygb", "Ttc9b", "Ipcef1"
-    )
-)
-
-makeSignMatrixPAGE(
-    sign_names = c("cell_type1", "cell_type2", "cell_type3"),
-    sign_list = sign_list
-)
-}
-\seealso{
-\code{\link{PAGEEnrich}}
-}
diff --git a/man/rankEnrich.Rd b/man/rankEnrich.Rd
deleted file mode 100644
index 7835c0802..000000000
--- a/man/rankEnrich.Rd
+++ /dev/null
@@ -1,39 +0,0 @@
-% Generated by roxygen2: do not edit by hand
-% Please edit documentation in R/spatial_enrichment.R
-\name{rankEnrich}
-\alias{rankEnrich}
-\title{rankEnrich}
-\usage{
-rankEnrich(...)
-}
-\arguments{
-\item{...}{
-  Arguments passed on to \code{\link[=runRankEnrich]{runRankEnrich}}
-  \describe{
-    \item{\code{gobject}}{Giotto object}
-    \item{\code{spat_unit}}{spatial unit}
-    \item{\code{feat_type}}{feature type}
-    \item{\code{sign_matrix}}{Matrix of signature genes for each cell type / process}
-    \item{\code{expression_values}}{expression values to use}
-    \item{\code{reverse_log_scale}}{reverse expression values from log scale}
-    \item{\code{logbase}}{log base to use if reverse_log_scale = TRUE}
-    \item{\code{output_enrichment}}{how to return enrichment output}
-    \item{\code{ties_method}}{how to handle rank ties}
-    \item{\code{p_value}}{calculate p-values (boolean, default = FALSE)}
-    \item{\code{n_times}}{number of permutations to calculate for p_value}
-    \item{\code{rbp_p}}{fractional binarization threshold (default = 0.99)}
-    \item{\code{num_agg}}{number of top genes to aggregate (default = 100)}
-    \item{\code{name}}{to give to spatial enrichment results, default = rank}
-    \item{\code{return_gobject}}{return giotto object}
-  }}
-}
-\value{
-gene enrichment
-}
-\description{
-Function to calculate gene signature enrichment scores per
-spatial position using a rank based approach.
-}
-\seealso{
-\code{\link{runRankEnrich}}
-}
diff --git a/man/runPAGEEnrich.Rd b/man/runPAGEEnrich.Rd
deleted file mode 100644
index 8d1ec1818..000000000
--- a/man/runPAGEEnrich.Rd
+++ /dev/null
@@ -1,100 +0,0 @@
-% Generated by roxygen2: do not edit by hand
-% Please edit documentation in R/spatial_enrichment.R
-\name{runPAGEEnrich}
-\alias{runPAGEEnrich}
-\title{runPAGEEnrich}
-\usage{
-runPAGEEnrich(
-  gobject,
-  spat_unit = NULL,
-  feat_type = NULL,
-  sign_matrix,
-  expression_values = c("normalized", "scaled", "custom"),
-  min_overlap_genes = 5,
-  reverse_log_scale = TRUE,
-  logbase = 2,
-  output_enrichment = c("original", "zscore"),
-  p_value = FALSE,
-  include_depletion = FALSE,
-  n_times = 1000,
-  max_block = 2e+07,
-  name = NULL,
-  verbose = TRUE,
-  return_gobject = TRUE
-)
-}
-\arguments{
-\item{gobject}{Giotto object}
-
-\item{spat_unit}{spatial unit}
-
-\item{feat_type}{feature type}
-
-\item{sign_matrix}{Matrix of signature genes for each cell type / process}
-
-\item{expression_values}{expression values to use}
-
-\item{min_overlap_genes}{minimum number of overlapping genes in sign_matrix
-required to calculate enrichment}
-
-\item{reverse_log_scale}{reverse expression values from log scale}
-
-\item{logbase}{log base to use if reverse_log_scale = TRUE}
-
-\item{output_enrichment}{how to return enrichment output}
-
-\item{p_value}{calculate p-values (boolean, default = FALSE)}
-
-\item{include_depletion}{calculate both enrichment and depletion}
-
-\item{n_times}{number of permutations to calculate for p_value}
-
-\item{max_block}{number of lines to process together (default = 20e6)}
-
-\item{name}{to give to spatial enrichment results, default = PAGE}
-
-\item{verbose}{be verbose}
-
-\item{return_gobject}{return giotto object}
-}
-\value{
-data.table with enrichment results
-}
-\description{
-Function to calculate gene signature enrichment scores per
-spatial position using PAGE.
-}
-\details{
-sign_matrix: a binary matrix with genes as row names and cell-types as
-column names.
-Alternatively a list of signature genes can be provided to
-makeSignMatrixPAGE, which will create the matrix for you. \cr
-
-The enrichment Z score is calculated by using method (PAGE) from
-Kim SY et al., BMC bioinformatics, 2005 as
-\eqn{Z = ((Sm – mu)*m^(1/2)) / delta}.
-For each gene in each spot, mu is the fold change values versus the mean
-expression and delta is the standard deviation. Sm is the mean fold change
-value of a specific marker gene set and  m is the size of a given marker
-gene set.
-}
-\examples{
-g <- GiottoData::loadGiottoMini("visium")
-sign_gene <- c(
-    "Bcl11b", "Lmo1", "F3", "Cnih3", "Ppp1r3c", "Rims2", "Gfap",
-    "Gjc3", "Chrna4", "Prkcd", "Prr18", "Grb14", "Tprn", "Clic1", "Olig2",
-    "Hrh3", "Tmbim1", "Carhsp1", "Tmem88b", "Ugt8a", "Arpp19", "Lamp5",
-    "Galnt6", "Hlf", "Hs3st2", "Tbr1", "Myl4", "Cygb", "Ttc9b", "Ipcef1"
-)
-
-sign_matrix <- matrix(rnorm(length(sign_gene) * 3, mean = 10),
-    nrow = length(sign_gene)
-)
-rownames(sign_matrix) <- sign_gene
-colnames(sign_matrix) <- c("cell_type1", "cell_type2", "cell_type3")
-
-runPAGEEnrich(gobject = g, sign_matrix = sign_matrix)
-}
-\seealso{
-\code{\link{makeSignMatrixPAGE}}
-}
diff --git a/man/runPAGEEnrich_OLD.Rd b/man/runPAGEEnrich_OLD.Rd
deleted file mode 100644
index e6e826493..000000000
--- a/man/runPAGEEnrich_OLD.Rd
+++ /dev/null
@@ -1,64 +0,0 @@
-% Generated by roxygen2: do not edit by hand
-% Please edit documentation in R/spatial_enrichment.R
-\name{runPAGEEnrich_OLD}
-\alias{runPAGEEnrich_OLD}
-\title{runPAGEEnrich_OLD}
-\usage{
-runPAGEEnrich_OLD(
-  gobject,
-  sign_matrix,
-  expression_values = c("normalized", "scaled", "custom"),
-  reverse_log_scale = TRUE,
-  logbase = 2,
-  output_enrichment = c("original", "zscore"),
-  p_value = FALSE,
-  n_times = 1000,
-  name = NULL,
-  return_gobject = TRUE
-)
-}
-\arguments{
-\item{gobject}{Giotto object}
-
-\item{sign_matrix}{Matrix of signature genes for each cell type / process}
-
-\item{expression_values}{expression values to use}
-
-\item{reverse_log_scale}{reverse expression values from log scale}
-
-\item{logbase}{log base to use if reverse_log_scale = TRUE}
-
-\item{output_enrichment}{how to return enrichment output}
-
-\item{p_value}{calculate p-values (boolean, default = FALSE)}
-
-\item{n_times}{number of permutations to calculate for p_value}
-
-\item{name}{to give to spatial enrichment results, default = PAGE}
-
-\item{return_gobject}{return giotto object}
-}
-\value{
-data.table with enrichment results
-}
-\description{
-Function to calculate gene signature enrichment scores per
-spatial position using PAGE.
-}
-\details{
-sign_matrix: a binary matrix with genes as row names and cell-types as
-column names.
-Alternatively a list of signature genes can be provided to
-makeSignMatrixPAGE, which will create the matrix for you. \cr
-
-The enrichment Z score is calculated by using method (PAGE) from
-Kim SY et al., BMC bioinformatics, 2005 as
-\eqn{Z = ((Sm – mu)*m^(1/2)) / delta}.
-For each gene in each spot, mu is the fold change values versus the mean
-expression and delta is the standard deviation. Sm is the mean fold change
-value of a specific marker gene set and  m is the size of a given marker
-gene set.
-}
-\seealso{
-\code{\link{makeSignMatrixPAGE}}
-}
diff --git a/man/screePlot.Rd b/man/screePlot.Rd
index 2f65fb4f9..a8477d3b1 100644
--- a/man/screePlot.Rd
+++ b/man/screePlot.Rd
@@ -8,7 +8,8 @@ screePlot(
   gobject,
   spat_unit = NULL,
   feat_type = NULL,
-  name = NULL,
+  dim_reduction_name = NULL,
+  name = deprecated(),
   expression_values = c("normalized", "scaled", "custom"),
   reduction = c("cells", "feats"),
   method = c("irlba", "exact", "random", "factominer"),
@@ -34,7 +35,9 @@ screePlot(
 
 \item{feat_type}{feature type (e.g. "rna", "dna", "protein")}
 
-\item{name}{name of PCA object if available}
+\item{dim_reduction_name}{name of PCA}
+
+\item{name}{deprecated}
 
 \item{expression_values}{expression values to use}
 
diff --git a/man/spatCellCellcomSpots.Rd b/man/spatCellCellcomSpots.Rd
index ac6e46389..dee251301 100644
--- a/man/spatCellCellcomSpots.Rd
+++ b/man/spatCellCellcomSpots.Rd
@@ -10,6 +10,7 @@ spatCellCellcomSpots(
   feat_type = NULL,
   ave_celltype_exp,
   spatial_network_name = "Delaunay_network",
+  spat_enr_name = "DWLS",
   cluster_column = "cell_ID",
   random_iter = 1000,
   feature_set_1,
@@ -39,6 +40,9 @@ spatCellCellcomSpots(
 \item{spatial_network_name}{spatial network to use for identifying
 interacting cells}
 
+\item{spat_enr_name}{name of spatial enrichment containing DWLS results.
+Default = \code{"DWLS"}}
+
 \item{cluster_column}{cluster column with cell type information}
 
 \item{random_iter}{number of iterations}
@@ -84,31 +88,33 @@ Statistical framework to identify if pairs of features
 expected based on a reshuffled null distribution of feature expression
 values in cells that are spatially in proximity to each other.
 \itemize{
- * LR_comb:Pair of ligand and receptor
- * lig_cell_type: cell type to assess expression level of ligand
- * lig_expr: average expression residual(observed - DWLS_predicted) of
- ligand in lig_cell_type
- * ligand: ligand name
- * rec_cell_type: cell type to assess expression level of receptor
- * rec_expr: average expression residual(observed - DWLS_predicted) of
- receptor in rec_cell_type
- * receptor: receptor name
- * LR_expr: combined average ligand and receptor expression residual
- * lig_nr: total number of cells from lig_cell_type that spatially interact
- with cells from rec_cell_type
- * rec_nr: total number of cells from rec_cell_type that spatially interact
- with cells from lig_cell_type
- * rand_expr: average combined ligand and receptor expression residual from
- random spatial permutations
- * av_diff: average difference between LR_expr and rand_expr over all random
- spatial permutations
- * sd_diff: (optional) standard deviation of the difference between LR_expr
- and rand_expr over all random spatial permutations
- * z_score: (optional) z-score
- * log2fc: LR_expr - rand_expr
- * pvalue: p-value
- * LR_cell_comb: cell type pair combination
- * p.adj: adjusted p-value
- * PI: significanc score: log2fc \* -log10(p.adj)
+\itemize{
+\item LR_comb:Pair of ligand and receptor
+\item lig_cell_type: cell type to assess expression level of ligand
+\item lig_expr: average expression residual(observed - DWLS_predicted) of
+ligand in lig_cell_type
+\item ligand: ligand name
+\item rec_cell_type: cell type to assess expression level of receptor
+\item rec_expr: average expression residual(observed - DWLS_predicted) of
+receptor in rec_cell_type
+\item receptor: receptor name
+\item LR_expr: combined average ligand and receptor expression residual
+\item lig_nr: total number of cells from lig_cell_type that spatially interact
+with cells from rec_cell_type
+\item rec_nr: total number of cells from rec_cell_type that spatially interact
+with cells from lig_cell_type
+\item rand_expr: average combined ligand and receptor expression residual from
+random spatial permutations
+\item av_diff: average difference between LR_expr and rand_expr over all random
+spatial permutations
+\item sd_diff: (optional) standard deviation of the difference between LR_expr
+and rand_expr over all random spatial permutations
+\item z_score: (optional) z-score
+\item log2fc: LR_expr - rand_expr
+\item pvalue: p-value
+\item LR_cell_comb: cell type pair combination
+\item p.adj: adjusted p-value
+\item PI: significanc score: log2fc \* -log10(p.adj)
+}
 }
 }