navinlabcode · WandeRum · Aug 23, 2026 · Aug 23, 2026
diff --git a/.gitignore b/.gitignore
@@ -0,0 +1,17 @@
+# Large test datasets (downloaded locally, not versioned)
+data/
+data2/
+data3/
+
+# Local test caches / logs / generated outputs
+tests_local/*.rds
+tests_local/*.log
+untitled_consensus_cluster/
+
+# R / RStudio
+.Rproj.user/
+.Rhistory
+.RData
+.Ruserdata
+*.o
+*.so
diff --git a/DESCRIPTION b/DESCRIPTION
@@ -1,6 +1,6 @@
 Package: CellTrek
 Title: Spatial Charting of Single Cell Transcriptomes in Tissues
-Version: 0.0.94
+Version: 0.0.95
 Authors@R: c(
     person("Runmin", 'Wei', email='wander1021@gmail.com', role=c('aut', 'cre')),
     person("Siyuan", 'He', email='SHe2@mdanderson.org', role=c('aut', 'ctb')) )
@@ -10,8 +10,8 @@ License: `use_mit_license()`, `use_gpl3_license()` or friends to pick a
 Encoding: UTF-8
 Roxygen: list(markdown = TRUE)
 RoxygenNote: 7.1.2
-Imports: 
-    akima,
+Imports:
+    interp,
     data.table,
     dbscan,
     dplyr,

diff --git a/NAMESPACE b/NAMESPACE
@@ -37,6 +37,6 @@ import(shiny)
 import(tibble)
 import(tidyr)
 import(visNetwork)
-importFrom(akima,interpp)
+importFrom(interp,interpp)
 importFrom(geometry,delaunayn)
 importFrom(packcircles,circleRepelLayout)
diff --git a/NEWS.md b/NEWS.md
@@ -0,0 +1,34 @@
+# CellTrek 0.0.95
+
+## Seurat v5 compatibility
+
+CellTrek now runs on Seurat v5 (tested with Seurat 5.3.0 / SeuratObject 5.2.0 on R 4.3).
+Verified end-to-end on VisiumV1 (brain_st_cortex + brain_sc) and VisiumV2 (stxBrain +
+allen_cortex) data through traint, celltrek, scoloc (KL/DT/KD), scoexp (cc) and run_kdist.
+
+* Replaced direct assay-slot access (`obj[[assay]]@data` / `@counts` / `@scale.data`)
+  with the version-agnostic `GetAssayData()` / `SetAssayData()` accessors throughout
+  `celltrek.R`, `scoexp.R` and `FindCorMarkers.R`. The old v3 idiom of moving counts into
+  the data slot now uses `SetAssayData(..., layer = 'data')`.
+* Made dimensional-reduction keys v5-compliant (`celltrekraw_`, `celltrek_`, `pca_`,
+  `umap_`, `tsne_`).
+* `CreateSeuratObject()` inputs are coerced with `as.matrix()` where needed.
+
+## Spatial coordinates
+
+* `traint()` resolves ST coordinates with the priority: `coord_df` > the image
+  `@coordinates[,coord_xy]` slot (VisiumV1, kept at original full-resolution scale) >
+  `GetTissueCoordinates()` (VisiumV2/FOV `x`/`y`). Preferring `@coordinates` avoids the
+  coordinate down-scaling that `GetTissueCoordinates()` applies to VisiumV1 objects
+  (which shrank spot distances and made the fixed `repel_r` over-disperse the charting).
+  Writes to `image@coordinates` are guarded to VisiumV1 (VisiumV2 has no such slot).
+* New `traint(coord_df = )` argument: supply ST spatial coordinates directly to bypass
+  differences/naming issues across image structures. Resolution priority is
+  `coord_df` > auto-detect > a clear error asking for `coord_df`. Rownames must cover all
+  ST cells/spots; columns named in `coord_xy` are used if present, otherwise the first two
+  columns are taken as (coord_x, coord_y).
+
+## Dependencies
+
+* Interpolation now uses `interp::interpp` instead of `akima::interpp` (akima is not
+  available on all platforms, e.g. macOS arm64). Updated DESCRIPTION and NAMESPACE.
diff --git a/R/FindCorMarkers.R b/R/FindCorMarkers.R
@@ -12,7 +12,7 @@
 #'
 #' @examples CorMarkers <- FindCorMarkers(srt_inp=test_seurat, assay='RNA', features=c('k_dist', 'sig_score'), method='spearman')
 FindCorMarkers <- function(srt_inp, assay='RNA', features, method='spearman') {
-  exp_dat <- as.matrix(srt_inp[[assay]]@data)
+  exp_dat <- as.matrix(GetAssayData(srt_inp, assay=assay, layer='data'))
   output_df <- data.frame(p_val=numeric(), cor=numeric(), p_val_adj=numeric(), feature=character(), gene=character())
   if (all(features %in% colnames(srt_inp@meta.data)))
   for (i in 1:length(features)) {

diff --git a/R/celltrek.R b/R/celltrek.R
diff --git a/R/scoexp.R b/R/scoexp.R
@@ -225,18 +225,18 @@ scoexp <- function(celltrek_inp, sigm=NULL, assay='RNA', gene_select=NULL, zero_
   if (is.null(sigm)) sigm <- celltrek_inp@images[[1]]@scale.factors$spot_dis
   if (is.null(gene_select)) {
     cat('gene filtering...\n')
-    feature_nz <- apply(celltrek_inp[[assay]]@data, 1, function(x) mean(x!=0)*100)
+    feature_nz <- apply(GetAssayData(celltrek_inp, assay=assay, layer='data'), 1, function(x) mean(x!=0)*100)
     features <- names(feature_nz)[feature_nz > zero_cutoff]
     cat(length(features), 'features after filtering...\n')
   } else if (length(gene_select) > 1) {
-    features <- intersect(gene_select, rownames(celltrek_inp[[assay]]@data))
+    features <- intersect(gene_select, rownames(GetAssayData(celltrek_inp, assay=assay, layer='data')))
     if (length(features)==0) stop('No genes in gene_select detected')
   }
   celltrek_inp <- Seurat::ScaleData(celltrek_inp, features=features)
   res <- list(gs=c(), cc=c(), rbfk=c(), wcor=c())
   dist_mat <- dist(celltrek_inp@meta.data[, c('coord_x', 'coord_y')]) %>% as.matrix
   kern_mat <- rbfk(dist_mat, sigm=sigm, zero_diag=F)
-  expr_mat <- t(as.matrix(celltrek_inp[[assay]]@scale.data))
+  expr_mat <- t(as.matrix(GetAssayData(celltrek_inp, assay=assay, layer='scale.data')))
   cat('Calculating spatial-weighted cross-correlation...\n')
   wcor_mat <- wcor(X=expr_mat, W=kern_mat, method=cor_method)
 

diff --git a/README.md b/README.md
@@ -8,6 +8,8 @@ In this tutorial, we will demonstrate the cell charting workflow based on the mo
 library(devtools)
 install_github("navinlabcode/CellTrek")
 ```
+
+> **Seurat v5 compatibility note.** This version has been updated to run on Seurat v5 (tested with Seurat 5.3.0 / SeuratObject 5.2.0 on R 4.3). Internally, direct assay-slot access (`obj[[assay]]@data`/`@counts`/`@scale.data`) was replaced with the version-agnostic `GetAssayData()`/`SetAssayData()` accessors, and dimensional-reduction keys were made v5-compliant. The interpolation step now depends on the `interp` package (a drop-in replacement for `akima`, which is not available on all platforms). The tutorial code below has likewise been updated for v5 (e.g. `HVFInfo()` instead of `@meta.features`).
 ## 2. Loading the packages and datasets (scRNA-seq and ST data)
 We start by loading the packages needed for the analyses. Please install them if you haven't.
 ``` r
@@ -124,20 +126,19 @@ Based on the CellTrek result, we can further investigate the co-expression patte
 L5 IT cells first are extracted from the charting result.
 ``` r
 brain_celltrek_l5 <- subset(brain_celltrek, subset=cell_type=='L5 IT')
-brain_celltrek_l5@assays$RNA@scale.data <- matrix(NA, 1, 1)
 brain_celltrek_l5$cluster <- gsub('L5 IT VISp ', '', brain_celltrek_l5$cluster)
 DimPlot(brain_celltrek_l5, group.by = 'cluster')
 ```
 ![](vignette_files/F6_L5IT.png)
 We select top 2000 variable genes (exclude mitochondrial, ribosomal and high-zero genes)
 ``` r
 brain_celltrek_l5 <- FindVariableFeatures(brain_celltrek_l5)
-vst_df <- brain_celltrek_l5@assays$RNA@meta.features %>% data.frame %>% mutate(id=rownames(.))
-nz_test <- apply(as.matrix(brain_celltrek_l5[['RNA']]@data), 1, function(x) mean(x!=0)*100)
+vst_df <- HVFInfo(brain_celltrek_l5) %>% data.frame %>% mutate(id=rownames(.))
+nz_test <- apply(as.matrix(GetAssayData(brain_celltrek_l5, assay='RNA', layer='data')), 1, function(x) mean(x!=0)*100)
 hz_gene <- names(nz_test)[nz_test<20]
 mt_gene <- grep('^Mt-', rownames(brain_celltrek_l5), value=T)
 rp_gene <- grep('^Rpl|^Rps', rownames(brain_celltrek_l5), value=T)
-vst_df <- vst_df %>% dplyr::filter(!(id %in% c(mt_gene, rp_gene, hz_gene))) %>% arrange(., -vst.variance.standardized)
+vst_df <- vst_df %>% dplyr::filter(!(id %in% c(mt_gene, rp_gene, hz_gene))) %>% arrange(., -variance.standardized)
 feature_temp <- vst_df$id[1:2000]
 ```
 We use scoexp to do the spatial-weighted gene co-expression analysis.

diff --git a/man/traint.Rd b/man/traint.Rd
diff --git a/tests_local/figs/1_celltrek_spatial.png b/tests_local/figs/1_celltrek_spatial.png
diff --git a/tests_local/figs/2_scoloc_KL_mst.png b/tests_local/figs/2_scoloc_KL_mst.png
diff --git a/tests_local/figs/3_scoexp_heatmap.png b/tests_local/figs/3_scoexp_heatmap.png
diff --git a/tests_local/figs/4_dataset2_spatial.png b/tests_local/figs/4_dataset2_spatial.png
diff --git a/tests_local/make_figs.R b/tests_local/make_figs.R
@@ -0,0 +1,60 @@
+suppressPackageStartupMessages({
+  library(CellTrek); library(Seurat); library(dplyr)
+  library(ggplot2); library(pheatmap); library(viridis)
+})
+set.seed(1)
+dir.create('tests_local/figs', showWarnings = FALSE)
+ct <- readRDS('tests_local/brain_celltrek_cache.rds')
+
+## ---- Fig 1: spatial charting (charted single cells on tissue coords) ----
+df <- ct@meta.data
+df$cell_type <- factor(df$cell_type, levels = sort(unique(df$cell_type)))
+p1 <- ggplot(df, aes(coord_y, max(coord_x)+min(coord_x)-coord_x, color = cell_type)) +
+  geom_point(size = 0.7, alpha = 0.85) +
+  coord_fixed() + theme_void() +
+  guides(color = guide_legend(override.aes = list(size = 3))) +
+  labs(title = paste0('CellTrek charting: ', nrow(df), ' single cells mapped to tissue'),
+       color = 'cell type')
+ggsave('tests_local/figs/1_celltrek_spatial.png', p1, width = 9, height = 7, dpi = 130)
+cat('fig1 done\n')
+
+## ---- Fig 2: SColoc (KL) consensus co-localization matrix ----
+glut <- c('L2/3 IT','L4','L5 IT','L5 PT','NP','L6 IT','L6 CT','L6b')
+g <- subset(ct, subset = cell_type %in% glut)
+sg <- CellTrek::scoloc(g, col_cell = 'cell_type', use_method = 'KL', eps = 1e-50, boot_n = 20)
+mst <- as.matrix(sg$mst_cons)
+labs <- glut; names(labs) <- make.names(glut)
+rownames(mst) <- colnames(mst) <- labs[rownames(mst)]
+mst_sym <- mst + t(mst)
+png('tests_local/figs/2_scoloc_KL_mst.png', width = 900, height = 800, res = 130)
+pheatmap(mst_sym, cluster_rows = TRUE, cluster_cols = TRUE,
+         color = viridis(100), display_numbers = TRUE, number_format = '%.2f',
+         main = 'SColoc (KL): MST consensus co-localization\n(glutamatergic neurons, 20 bootstraps)')
+dev.off()
+cat('fig2 done\n')
+
+## ---- Fig 3: SCoexp (cc) spatial co-expression heatmap on L5 IT ----
+l5 <- subset(ct, subset = cell_type == 'L5 IT')
+l5 <- FindVariableFeatures(l5, verbose = FALSE)
+vst_df <- HVFInfo(l5) %>% data.frame %>% mutate(id = rownames(.))
+nz <- apply(as.matrix(GetAssayData(l5, assay='RNA', layer='data')), 1, function(x) mean(x!=0)*100)
+hz <- names(nz)[nz < 20]
+mt <- grep('^Mt-', rownames(l5), value = TRUE)
+rp <- grep('^Rpl|^Rps', rownames(l5), value = TRUE)
+vst_df <- vst_df %>% filter(!(id %in% c(mt, rp, hz))) %>% arrange(-variance.standardized)
+feats <- head(na.omit(vst_df$id), 300)
+res <- CellTrek::scoexp(celltrek_inp = l5, assay = 'RNA', approach = 'cc',
+                        gene_select = feats, sigm = 140, avg_cor_min = .3,
+                        zero_cutoff = 3, min_gen = 15, max_gen = 200, maxK = 6, k = 6, reps = 20)
+k_df <- do.call(rbind, lapply(seq_along(res$gs), function(i)
+  data.frame(gene = res$gs[[i]], Module = paste0('K', i))))
+rownames(k_df) <- k_df$gene; k_df$gene <- NULL
+gord <- rownames(k_df)
+png('tests_local/figs/3_scoexp_heatmap.png', width = 850, height = 800, res = 130)
+pheatmap(res$wcor[gord, gord], clustering_method = 'ward.D2',
+         annotation_row = k_df, show_rownames = FALSE, show_colnames = FALSE,
+         treeheight_row = 12, treeheight_col = 12, fontsize = 8,
+         color = viridis(10), main = 'SCoexp: L5 IT spatial co-expression modules')
+dev.off()
+cat('fig3 done | modules:', length(res$gs), '| sizes:', paste(sapply(res$gs, length), collapse=','), '\n')
+cat('ALL FIGS DONE\n')
diff --git a/tests_local/regress_coords.R b/tests_local/regress_coords.R
@@ -0,0 +1,28 @@
+## Regression check after the coordinate-scale fix: traint must resolve coords for
+## both VisiumV1 (brain_st_cortex, via @coordinates) and VisiumV2 (stxBrain, via GTC x/y).
+suppressPackageStartupMessages({library(CellTrek); library(Seurat); library(dplyr)})
+set.seed(1)
+mini_sc <- function(sc, per=25){ md<-sc@meta.data; k<-unlist(lapply(split(rownames(md), md$cell_type), function(x) sample(x, min(length(x),per)))); sc[,k] }
+
+## ---- VisiumV1 ----
+st1 <- readRDS('data/brain_st_cortex.rds'); sc1 <- readRDS('data/brain_sc.rds')
+st1 <- RenameCells(st1, new.names=make.names(Cells(st1))); sc1 <- RenameCells(sc1, new.names=make.names(Cells(sc1)))
+st1 <- st1[, sample(Cells(st1), 300)]; sc1 <- mini_sc(sc1, 25)
+t1 <- CellTrek::traint(st1, sc1, st_assay='Spatial', sc_assay='RNA', cell_names='cell_type')
+c1 <- na.omit(t1@meta.data[t1$type=='st', c('coord_x','coord_y')])
+cat('VisiumV1: image', class(st1@images[[1]]), '| coord range', paste(round(range(c1$coord_x)),collapse='-'),
+    '| spot_dis', round(median(dbscan::kNN(c1,k=6)$dist)), '| st coords set', all(!is.na(t1$coord_x[t1$type=='st'])), '\n')
+
+## ---- VisiumV2 ----
+e <- new.env(); load('data2/anterior1.rda', envir=e); st2 <- UpdateSeuratObject(get('anterior1', envir=e))
+st2 <- NormalizeData(st2, verbose=FALSE) %>% FindVariableFeatures(nfeatures=2000, verbose=FALSE)
+st2 <- RenameCells(st2, new.names=make.names(Cells(st2)))
+sc2 <- readRDS('data2/allen_cortex.rds'); sc2$cell_type <- as.character(sc2$subclass)
+sc2 <- mini_sc(sc2, 25)
+sc2 <- NormalizeData(sc2,verbose=FALSE) %>% FindVariableFeatures(nfeatures=2000,verbose=FALSE) %>% ScaleData(verbose=FALSE) %>% RunPCA(npcs=30,verbose=FALSE) %>% RunUMAP(dims=1:30,verbose=FALSE)
+sc2 <- RenameCells(sc2, new.names=make.names(Cells(sc2)))
+t2 <- CellTrek::traint(st2, sc2, st_assay='Spatial', sc_assay='RNA', cell_names='cell_type')
+c2 <- na.omit(t2@meta.data[t2$type=='st', c('coord_x','coord_y')])
+cat('VisiumV2: image', class(st2@images[[1]]), '| coord range', paste(round(range(c2$coord_x)),collapse='-'),
+    '| spot_dis', round(median(dbscan::kNN(c2,k=6)$dist)), '| st coords set', all(!is.na(t2$coord_x[t2$type=='st'])), '\n')
+cat('== REGRESS OK ==\n')
diff --git a/tests_local/run_dataset2.R b/tests_local/run_dataset2.R
@@ -0,0 +1,69 @@
+## Second paired dataset test on Seurat v5:
+##   ST  = stxBrain anterior1 (10x Visium, updates to Assay5 + VisiumV2 image)
+##   scRNA = allen_cortex reference (subclass labels)
+## This exercises the VisiumV2 (v5 image class) code path, which the first
+## dataset (brain_st_cortex, VisiumV1) did not.
+suppressPackageStartupMessages({
+  library(CellTrek); library(Seurat); library(dplyr); library(magrittr)
+  library(ggplot2)
+})
+set.seed(1)
+cat('Seurat', as.character(packageVersion('Seurat')), '\n')
+dir.create('tests_local/figs', showWarnings = FALSE)
+
+## ---- ST: stxBrain anterior1 ----
+e <- new.env(); load('data2/anterior1.rda', envir = e)
+st <- UpdateSeuratObject(get('anterior1', envir = e))
+cat('ST anterior1:', ncol(st), 'spots | image class:', class(st@images[[1]]), '\n')
+st <- NormalizeData(st, verbose = FALSE) %>% FindVariableFeatures(nfeatures = 2000, verbose = FALSE)
+st <- RenameCells(st, new.names = make.names(Cells(st)))
+
+## ---- scRNA: allen_cortex, subsample + preprocess ----
+sc <- readRDS('data2/allen_cortex.rds')
+sc$cell_type <- as.character(sc$subclass)
+set.seed(42)
+md <- sc@meta.data
+keep <- unlist(lapply(split(rownames(md), md$cell_type), function(x) sample(x, min(length(x), 70))))
+sc <- subset(sc, cells = keep)
+sc <- NormalizeData(sc, verbose = FALSE) %>%
+  FindVariableFeatures(nfeatures = 2000, verbose = FALSE) %>%
+  ScaleData(verbose = FALSE) %>% RunPCA(npcs = 30, verbose = FALSE) %>%
+  RunUMAP(dims = 1:30, verbose = FALSE)
+sc <- RenameCells(sc, new.names = make.names(Cells(sc)))
+cat('scRNA allen:', ncol(sc), 'cells |', length(unique(sc$cell_type)), 'cell types\n')
+
+## ---- Step 1: traint ----
+cat('\n==== traint (VisiumV2 ST) ====\n')
+tr <- CellTrek::traint(st_data = st, sc_data = sc, st_assay = 'Spatial',
+                       sc_assay = 'RNA', cell_names = 'cell_type')
+cat('traint:', ncol(tr), 'cells | coord_x set:',
+    all(!is.na(tr$coord_x[tr$type == 'st'])), '\n')
+
+## ---- Step 2: celltrek ----
+cat('\n==== celltrek ====\n')
+res <- CellTrek::celltrek(st_sc_int = tr, int_assay = 'traint', sc_data = sc,
+                          sc_assay = 'RNA', reduction = 'pca', intp = TRUE,
+                          intp_pnt = 5000, intp_lin = FALSE, nPCs = 30, ntree = 1000,
+                          dist_thresh = 0.55, top_spot = 5, spot_n = 5,
+                          repel_r = 20, repel_iter = 20, keep_model = FALSE)$celltrek
+cat('celltrek:', ncol(res), 'charted cells | coords in meta:',
+    all(c('coord_x','coord_y') %in% colnames(res@meta.data)), '\n')
+saveRDS(res, 'tests_local/dataset2_celltrek_cache.rds')
+
+## ---- Fig: spatial charting ----
+df <- res@meta.data
+df$cell_type <- factor(df$cell_type, levels = sort(unique(df$cell_type)))
+p <- ggplot(df, aes(coord_y, max(coord_x)+min(coord_x)-coord_x, color = cell_type)) +
+  geom_point(size = 0.6, alpha = 0.85) + coord_fixed() + theme_void() +
+  guides(color = guide_legend(override.aes = list(size = 3))) +
+  labs(title = paste0('Dataset 2 (stxBrain anterior + allen_cortex): ',
+                      nrow(df), ' cells charted'), color = 'subclass')
+ggsave('tests_local/figs/4_dataset2_spatial.png', p, width = 9, height = 7, dpi = 130)
+cat('spatial fig saved\n')
+
+## ---- Step 3: scoloc (DT, no external dep) ----
+cat('\n==== scoloc (DT) ====\n')
+sg <- CellTrek::scoloc(res, col_cell = 'cell_type', use_method = 'DT', boot_n = 10)
+cat('scoloc mst_cons:', paste(dim(sg$mst_cons), collapse = 'x'), '\n')
+
+cat('\n==== DATASET 2 ALL STAGES COMPLETED ====\n')