navinlabcode / CellTrek

103 stars 18 forks source link

Error in rownames(st_data[[st_assay]]@data) : no slot of name "data" for this object of class "Assay5", Possible Seurat v5.0 error #39

Open tsoruu opened 9 months ago

tsoruu commented 9 months ago

Hiya, got a problem here

traint_result <- CellTrek::traint(st_data = ST,
                                  sc_data = SC,
                                  sc_assay = 'RNA',
                                  norm = 'SCT',
                                  cell_names = 'subclass')
Finding transfer anchors... 
No variable features found for object2 in the object.list. Running FindVariableFeatures ...
Calculating gene variances
0%   10   20   30   40   50   60   70   80   90   100%
[----|----|----|----|----|----|----|----|----|----|
**************************************************|
Calculating feature variances of standardized and clipped values
0%   10   20   30   40   50   60   70   80   90   100%
[----|----|----|----|----|----|----|----|----|----|
**************************************************|
Error in rownames(st_data[[st_assay]]@data) : 
  no slot of name "data" for this object of class "Assay5"

I think the error is caused by Seurat v5 object's format change. This is the format of my object, which uses Seurat v5 image The count and data are wrapped in layer object(?)/variable(?). I wonder if someone has a tool to downgrade Seurat object, since my preprocessing code didn't work on the latest v4...

p.s: Tried to just copy and paste the count and data outside the layer object(?)/variable(?), didn't work...

sessionInfo()

R version 4.3.2 (2023-10-31)
Platform: aarch64-unknown-linux-gnu (64-bit)
Running under: Ubuntu 22.04.3 LTS

Matrix products: default
BLAS:   /usr/lib/aarch64-linux-gnu/openblas-pthread/libblas.so.3 
LAPACK: /usr/lib/aarch64-linux-gnu/openblas-pthread/libopenblasp-r0.3.20.so;  LAPACK version 3.10.0

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

time zone: Asia/Jakarta
tzcode source: system (glibc)

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

other attached packages:
[1] ConsensusClusterPlus_1.66.0 viridis_0.6.4               viridisLite_0.4.2           Seurat_5.0.1               
[5] SeuratObject_5.0.1          sp_2.1-2                    dplyr_1.1.4                 CellTrek_0.0.94            

loaded via a namespace (and not attached):
  [1] RColorBrewer_1.1-3     rstudioapi_0.15.0      jsonlite_1.8.8         magrittr_2.0.3        
  [5] spatstat.utils_3.0-4   vctrs_0.6.5            ROCR_1.0-11            spatstat.explore_3.2-5
  [9] rstatix_0.7.2          htmltools_0.5.7        dynamicTreeCut_1.63-1  broom_1.0.5           
 [13] sctransform_0.4.1      parallelly_1.36.0      KernSmooth_2.23-22     htmlwidgets_1.6.4     
 [17] ica_1.0-3              plyr_1.8.9             plotly_4.10.3          zoo_1.8-12            
 [21] igraph_1.6.0           mime_0.12              lifecycle_1.0.4        pkgconfig_2.0.3       
 [25] Matrix_1.6-4           R6_2.5.1               fastmap_1.1.1          magic_1.6-1           
 [29] fitdistrplus_1.1-11    future_1.33.1          shiny_1.8.0            digest_0.6.33         
 [33] colorspace_2.1-0       patchwork_1.1.3        tensor_1.5             RSpectra_0.16-1       
 [37] irlba_2.3.5.1          akima_0.6-3.4          ggpubr_0.6.0           philentropy_0.8.0     
 [41] progressr_0.14.0       fansi_1.0.6            spatstat.sparse_3.0-3  httr_1.4.7            
 [45] polyclip_1.10-6        abind_1.4-5            compiler_4.3.2         backports_1.4.1       
 [49] carData_3.0-5          fastDummies_1.7.3      ggsignif_0.6.4         MASS_7.3-60           
 [53] tools_4.3.2            lmtest_0.9-40          httpuv_1.6.13          future.apply_1.11.1   
 [57] goftest_1.2-3          glue_1.6.2             dbscan_1.1-12          DiagrammeR_1.0.10     
 [61] nlme_3.1-163           promises_1.2.1         grid_4.3.2             Rtsne_0.17            
 [65] cluster_2.1.4          reshape2_1.4.4         generics_0.1.3         gtable_0.3.4          
 [69] spatstat.data_3.0-3    tidyr_1.3.0            data.table_1.14.10     car_3.1-2             
 [73] utf8_1.2.4             BiocGenerics_0.48.1    spatstat.geom_3.2-7    RcppAnnoy_0.0.21      
 [77] ggrepel_0.9.4          RANN_2.6.1             pillar_1.9.0           stringr_1.5.1         
 [81] spam_2.10-0            RcppHNSW_0.5.0         later_1.3.2            splines_4.3.2         
 [85] lattice_0.22-5         survival_3.5-7         deldir_2.0-2           tidyselect_1.2.0      
 [89] miniUI_0.1.1.1         pbapply_1.7-2          gridExtra_2.3          scattermore_1.2       
 [93] Biobase_2.62.0         matrixStats_1.2.0      visNetwork_2.1.2       stringi_1.8.3         
 [97] lazyeval_0.2.2         codetools_0.2-19       data.tree_1.1.0        tibble_3.2.1          
[101] packcircles_0.3.6      cli_3.6.2              uwot_0.1.16            geometry_0.4.7        
[105] xtable_1.8-4           reticulate_1.34.0      randomForestSRC_3.2.3  munsell_0.5.0         
[109] Rcpp_1.0.11            globals_0.16.2         spatstat.random_3.2-2  png_0.1-8             
[113] fastcluster_1.2.3      parallel_4.3.2         ellipsis_0.3.2         ggplot2_3.4.4         
[117] dotCall64_1.1-1        listenv_0.9.0          scales_1.3.0           ggridges_0.5.5        
[121] leiden_0.4.3.1         purrr_1.0.2            rlang_1.1.2            cowplot_1.1.2
tsoruu commented 8 months ago

Got some progress on this issue, I've tried converting the object to a similar format to v3 using https://github.com/satijalab/seurat/issues/7409#issuecomment-1584041598 code, but currently encounter this problem


> traint_result <- traint(st_data = ST,
+                         sc_data = SC,
+                         sc_assay = 'RNA',
+                         st_assay = 'Spatial',
+                         cell_names = 'subclass',
+                         norm = 'LogNormalize')
Finding transfer anchors... 
No variable features found for object2 in the object.list. Running FindVariableFeatures ...
Calculating gene variances
0%   10   20   30   40   50   60   70   80   90   100%
[----|----|----|----|----|----|----|----|----|----|
**************************************************|
Calculating feature variances of standardized and clipped values
0%   10   20   30   40   50   60   70   80   90   100%
[----|----|----|----|----|----|----|----|----|----|
**************************************************|
Using 2000 features for integration... 
Running CCA
Merging objects
Finding neighborhoods
Finding anchors
    Found 6443 anchors
Data transfering... 
Finding integration vectors
Finding integration vector weights
0%   10   20   30   40   50   60   70   80   90   100%
[----|----|----|----|----|----|----|----|----|----|
**************************************************|
Transfering 2000 features onto reference data
Creating new Seurat object... 
Warning: Data is of class data.frame. Coercing to dgCMatrix.
Error in traint(st_data = ST, sc_data = SC, sc_assay = "RNA", st_assay = "Spatial",  : 
  no slot of name "counts" for this object of class "Assay5"
In addition: Warning message:
Adding image data that isn't associated with any assays 
> traceback()
1: traint(st_data = ST, sc_data = SC, sc_assay = "RNA", st_assay = "Spatial", 
       cell_names = "subclass", norm = "LogNormalize") ```

While the Seurat spatial object that I've converted to v3/4-like format is something like this...
![image](https://github.com/navinlabcode/CellTrek/assets/83799574/c8a219f3-8d05-4f2f-a5bc-4116a33ae4d7)

The traceback is kind of weird there, but according to the error message, it seems that there isn't a count in class "Assay5", as you can see, there is no assay5 again, rather it's now just assay, indicating it's conforming with the v3/4 object. My guess would be after it tries to create a new Seurat object, it tries to check whether the object is successfully created or not by checking the count and data, but the new object conforms with the v5 format, making it unable to find the count. Is there any way to set the default Seurat object creation to conform with the v3/4 format?
Mapleljs commented 8 months ago

I am facing the same question as you. Have you known how to solve it?

xieaq commented 7 months ago

Same question :(

AndoneguiElguera commented 6 months ago

The issue is associated with object sctructure, I can´t change from Assay5 to Assay3 but I changed the functions code of CellTrek using the follow code:

#### CellTrek Functions ####
# traint ####
traint <- function (
    st_data, sc_data, st_assay='Spatial', sc_assay='scint', 
    norm='LogNormalize', nfeatures=2000, 
    cell_names='cell_names', coord_xy=c('imagerow', 'imagecol'), 
    gene_kept=NULL, ...) {
st_data$id <- names(st_data$orig.ident)
sc_data$id <- names(sc_data$orig.ident)
sc_data$cell_names <- make.names(sc_data@meta.data[, cell_names])
st_data$type <- 'st'
sc_data$type <- 'sc'
st_data$coord_x <- st_data@images[[1]]@coordinates[, coord_xy[1]]
st_data$coord_y <- st_data@images[[1]]@coordinates[, coord_xy[2]]
DefaultAssay(st_data) <- st_assay
DefaultAssay(sc_data) <- sc_assay

cat('Finding transfer anchors... \n')
st_idx <- st_data$id
sc_idx <- sc_data$id

## Integration features ##
sc_st_list <- list(st_data=st_data, sc_data=sc_data)
sc_st_features <- Seurat::SelectIntegrationFeatures(
  sc_st_list, nfeatures=nfeatures)
if (!is.null(gene_kept)) {
  sc_st_features <- union(sc_st_features, gene_kept)
}

sc_st_features <- 
  sc_st_features[(sc_st_features %in% rownames(st_data[[st_assay]]@data)) & 
                   (sc_st_features %in% rownames(sc_data[[sc_assay]]@data))]
cat('Using', length(sc_st_features), 'features for integration... \n')
###

sc_st_anchors <- Seurat::FindTransferAnchors(
  reference = sc_data, query = st_data, 
  reference.assay = sc_assay, query.assay = st_assay,
  normalization.method = norm, features = sc_st_features, 
  reduction = 'cca')

cat('Data transfering... \n')
st_data_trans <- Seurat::TransferData(
  anchorset = sc_st_anchors, 
  refdata = GetAssayData(sc_data)[sc_st_features, ], 
  weight.reduction = 'cca')

st_data@assays$transfer <- st_data_trans

cat('Creating new Seurat object... \n')
sc_st_meta <- dplyr::bind_rows(st_data@meta.data, sc_data@meta.data)
counts_temp <- cbind(
  data.frame(st_data[['transfer']]@data),
  data.frame(sc_data[[sc_assay]]@data[sc_st_features, ] %>% data.frame))
rownames(sc_st_meta) <- make.names(sc_st_meta$id)
colnames(counts_temp) <- make.names(sc_st_meta$id)
sc_st_int <- CreateSeuratObject(
  counts = counts_temp, 
  assay = 'traint', 
  meta.data = sc_st_meta, project = "SeuratProject"
)

cat('Scaling -> PCA -> UMAP... \n')
sc_st_int <- ScaleData(
  sc_st_int, 
  features = sc_st_features, 
  assay = 'traint', layer = 'counts') %>%
  RunPCA(features = sc_st_features)
sc_st_int <- RunUMAP(sc_st_int, dims = 1:30)
sc_st_int@images <- st_data@images
library(magrittr)
sc_st_int@images[[1]]@coordinates <- data.frame(
  imagerow=sc_st_int@meta.data$coord_x,
  imagecol=sc_st_int@meta.data$coord_y) %>% 
  set_rownames(rownames(sc_st_int@meta.data))

return (sc_st_int)
}

# celltrek ####
celltrek <- function (
    st_sc_int, int_assay='traint', sc_data=NULL, sc_assay='RNA', 
    reduction='pca', intp=T, intp_pnt=10000, intp_lin=F, nPCs=30, ntree=1000, 
    dist_thresh=.4, top_spot=10, spot_n=10, repel_r=5, repel_iter=10, 
    keep_model=F, ...) {

dist_res <- celltrek_dist(
  st_sc_int=st_sc_int, int_assay=int_assay, reduction=reduction, intp=intp, 
  intp_pnt=intp_pnt, intp_lin=intp_lin, nPCs=nPCs, ntree=ntree, keep_model=T)

spot_dis_intp <- mediaspot_dis_intp <- mediaspot_dis_intp <- median(
  unlist(dbscan::kNN(dist_res$coord_df[, c('coord_x', 'coord_y')], k=4)$dist))

if (is.null(repel_r)) {repel_r=spot_dis_intp/4}

sc_coord_list <- celltrek_chart(
  dist_mat=dist_res$celltrek_dist, coord_df=dist_res$coord_df, 
  dist_cut=ntree*dist_thresh, top_spot=top_spot, spot_n=spot_n, 
  repel_r=repel_r, repel_iter=repel_iter)

sc_coord_raw <- sc_coord_list[[1]]
sc_coord <- sc_coord_list[[2]]

cat('Creating Seurat Object... \n')
if (!is.null(sc_data)) {
  cat('sc data...')
  sc_data$id <- Seurat::Cells(sc_data)

  sc_out <- CreateSeuratObject(
    counts=sc_data[[sc_assay]]@data[, sc_coord$id_raw] %>% 
      set_colnames(sc_coord$id_new), 
    project='celltrek', 
    assay=sc_assay, 
    meta.data=sc_data@meta.data[sc_coord$id_raw, ] %>% 
      dplyr::rename(id_raw=id) %>% 
      mutate(id_new=sc_coord$id_new) %>% 
      set_rownames(sc_coord$id_new))

  sc_out@meta.data <- dplyr::left_join(
    sc_out@meta.data, sc_coord) %>% 
    data.frame %>% 
    set_rownames(sc_out$id_new)

  sc_coord_raw_df <- CreateDimReducObject(
    embeddings=sc_coord_raw %>% 
      dplyr::mutate(
        coord1=coord_y, coord2=max(coord_x)+min(coord_x)-coord_x
      ) %>%
      dplyr::select(c(coord1, coord2)) %>% 
      set_rownames(sc_coord_raw$id_new) %>% 
      as.matrix, assay=sc_assay, key='celltrek_raw')

  sc_coord_dr <- CreateDimReducObject(
    embeddings=sc_coord %>% 
      dplyr::mutate(
        coord1=coord_y, coord2=max(coord_x)+min(coord_x)-coord_x) %>% 
      dplyr::select(c(coord1, coord2)) %>% 
      set_rownames(sc_coord$id_new) %>% 
      as.matrix, assay=sc_assay, key='celltrek')

  sc_out@reductions$celltrek <- sc_coord_dr
  sc_out@reductions$celltrek_raw <- sc_coord_raw_df
  if ('pca' %in% names(sc_data@reductions)) {
    sc_pca_dr <- CreateDimReducObject(
      embeddings=sc_data@reductions$pca@cell.embeddings[sc_coord$id_raw, ] %>% 
        set_rownames(sc_coord$id_new) %>% as.matrix, assay=sc_assay, key='pca'
    )
    sc_out@reductions$pca <- sc_pca_dr
  }
  if ('umap' %in% names(sc_data@reductions)) {
    sc_umap_dr <- CreateDimReducObject(
      embeddings=sc_data@reductions$umap@cell.embeddings[sc_coord$id_raw,] %>% 
        set_rownames(sc_coord$id_new) %>% as.matrix, assay=sc_assay, 
      key='umap'
    )
    sc_out@reductions$umap <- sc_umap_dr
  }
  if ('tsne' %in% names(sc_data@reductions)) {
    sc_tsne_dr <- CreateDimReducObject(
      embeddings=sc_data@reductions$tsne@cell.embeddings[sc_coord$id_raw, ] %>% 
        set_rownames(sc_coord$id_new) %>% 
        as.matrix, assay=sc_assay, key='tsne')
    sc_out@reductions$tsne <- sc_tsne_dr
  }
} else {
  cat('no sc data...')
  sc_out <- CreateSeuratObject(
    counts=st_sc_int[[int_assay]]@data[, sc_coord$id_raw] %>% 
      set_colnames(sc_coord$id_new), 
    project='celltrek', assay=int_assay, 
    meta.data=st_sc_int@meta.data[sc_coord$id_raw, ] %>%
                                 dplyr::rename(id_raw=id) %>%
                                 mutate(id_new=sc_coord$id_new) %>%
                                 set_rownames(sc_coord$id_new)
    )
  sc_out$coord_x <- sc_coord$coord_x[match(sc_coord$id_new, sc_out$id_new)]
  sc_out$coord_y <- sc_coord$coord_y[match(sc_coord$id_new, sc_out$id_new)]

  sc_out[[int_assay]]@scale.data <- 
    st_sc_int[[int_assay]]@scale.data[, sc_coord$id_raw] %>% 
    set_colnames(sc_coord$id_new)
  sc_coord_raw_df <- CreateDimReducObject(
    embeddings=sc_coord_raw %>% 
      dplyr::mutate(
        coord1=coord_y, coord2=max(coord_x)+min(coord_x)-coord_x) %>%
      dplyr::select(c(coord1, coord2)) %>% 
      set_rownames(sc_coord_raw$id_new) %>% 
      as.matrix, assay=sc_assay, key='celltrek_raw'
    )
  sc_coord_dr <- CreateDimReducObject(
    embeddings = sc_coord %>% 
      dplyr::mutate(
        coord1=coord_y, coord2=max(coord_x)+min(coord_x)-coord_x) %>% 
      dplyr::select(c(coord1, coord2)) %>% 
      set_rownames(sc_coord$id_new) %>% 
      as.matrix, assay=int_assay, key='celltrek'
    )
  sc_out@reductions$celltrek <- sc_coord_dr
  sc_out@reductions$celltrek_raw <- sc_coord_raw_df
  if ('pca' %in% names(st_sc_int@reductions)) {
    sc_pca_dr <- CreateDimReducObject(
      embeddings=st_sc_int@reductions$pca@cell.embeddings[sc_coord$id_raw, ] %>% 
        set_rownames(sc_coord$id_new) %>%
        as.matrix, assay=int_assay, key='pca'
      )
    sc_out@reductions$pca <- sc_pca_dr
  }
  if ('umap' %in% names(st_sc_int@reductions)) {
    sc_umap_dr <- CreateDimReducObject(
      embeddings=st_sc_int@reductions$umap@cell.embeddings[sc_coord$id_raw,] %>% 
        set_rownames(sc_coord$id_new) %>%
        as.matrix, assay=int_assay, key='umap'
      )
    sc_out@reductions$umap <- sc_umap_dr
  }
  if ('tsne' %in% names(st_sc_int@reductions)) {
    sc_tsne_dr <- CreateDimReducObject(
      embeddings=st_sc_int@reductions$tsne@cell.embeddings[sc_coord$id_raw,] %>%
        set_rownames(sc_coord$id_new) %>%
        as.matrix, assay=int_assay, key='tsne'
      )
    sc_out@reductions$tsne <- sc_tsne_dr
  }
}
sc_out@images <- st_sc_int@images
sc_out@images[[1]]@assay <- DefaultAssay(sc_out)
sc_out@images[[1]]@coordinates <- data.frame(
  imagerow=sc_coord$coord_x, imagecol=sc_coord$coord_y) %>%
  set_rownames(sc_coord$id_new)
sc_out@images[[1]]@scale.factors$spot_dis <- dist_res$spot_d
sc_out@images[[1]]@scale.factors$spot_dis_intp <- spot_dis_intp

output <- list(celltrek=sc_out)
if (keep_model) {
  output[[length(output)+1]] <- dist_res$model
  names(output)[length(output)] <- 'model'
}
return(output)
}

# celltrek_dist ####
celltrek_dist <- function (
    st_sc_int, int_assay='traint', reduction='pca', intp = T, intp_pnt=10000, 
    intp_lin=F, nPCs=30, ntree=1000, keep_model=T) {
  DefaultAssay(st_sc_int) <- int_assay
  kNN_dist <- dbscan::kNN(
    na.omit(st_sc_int@meta.data[, c('coord_x', 'coord_y')]), k=6)$dist
  spot_dis <- median(kNN_dist) %>% round
  cat('Distance between spots is:', spot_dis, '\n')

  st_sc_int$id <- names(st_sc_int$orig.ident)
  st_idx <- st_sc_int$id[st_sc_int$type=='st']
  sc_idx <- st_sc_int$id[st_sc_int$type=='sc']
  meta_df <- data.frame(st_sc_int@meta.data)

  st_sc_int_pca <- 
    st_sc_int@reductions[[reduction]]@cell.embeddings[, 1:nPCs] %>% 
    data.frame %>%
    mutate(id=st_sc_int$id, type=st_sc_int$type, class=st_sc_int$cell_names,
           coord_x=st_sc_int$coord_x, coord_y=st_sc_int$coord_y)
  st_pca <- st_sc_int_pca %>% 
    dplyr::filter(type=='st') %>% 
    dplyr::select(-c(id:class))

  ## Interpolation ##
  ## Uniform sampling ##
  if (intp) {
    cat ('Interpolating...\n')
    spot_ratio <- intp_pnt/nrow(st_pca)
    st_intp_df <- apply(st_pca[, c('coord_x', 'coord_y')], 1, function(row_x) {
      runif_test <- runif(1)
      if (runif_test < spot_ratio%%1) {
        theta <- runif(ceiling(spot_ratio), 0, 2*pi)
        alpha <- sqrt(runif(ceiling(spot_ratio), 0, 1))
        coord_x <- row_x[1] + (spot_dis/2)*sin(theta)*alpha
        coord_y <- row_x[2] + (spot_dis/2)*cos(theta)*alpha
      } else {
        theta <- runif(floor(spot_ratio), 0, 2*pi)
        alpha <- sqrt(runif(floor(spot_ratio), 0, 1))
        coord_x <- row_x[1] + (spot_dis/2)*sin(theta)*alpha
        coord_y <- row_x[2] + (spot_dis/2)*cos(theta)*alpha
      }
      data.frame(coord_x, coord_y)
    }) %>% Reduce(rbind, .)

    st_intp_df <- apply(st_pca[, 1:nPCs], 2, function(col_x) {
      akima::interpp(x=st_pca$coord_x, y=st_pca$coord_y, z=col_x,
                     linear=intp_lin, xo=st_intp_df$coord_x, 
                     yo=st_intp_df$coord_y) %>%
        magrittr::extract2('z')
    }) %>% data.frame(., id='X', type='st_intp', st_intp_df) %>% na.omit
    st_intp_df$id <- make.names(st_intp_df$id, unique = T)
    st_sc_int_pca <- bind_rows(st_sc_int_pca, st_intp_df)
  }

  cat('Random Forest training... \n')
  ## Training on ST ##
  data_train <- st_sc_int_pca %>% 
    dplyr::filter(type=='st') %>% 
    dplyr::select(-c(id:class))
  rf_train <- randomForestSRC::rfsrc(
    Multivar(coord_x, coord_y) ~ ., data_train, block.size=5, ntree=ntree)

  cat('Random Forest prediction...  \n')
  ## Testing on all ##
  data_test <- st_sc_int_pca
  rf_pred <- randomForestSRC::predict.rfsrc(
    rf_train, newdata=data_test[, c(1:nPCs)], distance='all')

  cat('Making distance matrix... \n')
  rf_pred_dist <- 
    rf_pred$distance[data_test$type=='sc', data_test$type!='sc'] %>%
    set_rownames(data_test$id[data_test$type=='sc']) %>% 
    set_colnames(data_test$id[data_test$type!='sc'])

  output <- list()
  output$spot_d <- spot_dis
  output$celltrek_dist <- rf_pred_dist
  output$coord_df <- st_sc_int_pca[, c('id', 'type', 'coord_x', 'coord_y')] %>%
    dplyr::filter(type!='sc') %>% 
    magrittr::set_rownames(.$id) %>% 
    dplyr::select(-id)
  if (keep_model) {
    output$model <- rf_train
  }
  return (output)
}

# celltrek_chart ####
celltrek_chart <- function (dist_mat, coord_df, dist_cut=500, top_spot=10, 
                            spot_n=10, repel_r=5, repel_iter=10) {
  cat('Making graph... \n')
  dist_mat[dist_mat>dist_cut] <- NA
  dist_mat_dt <- data.table::data.table(Var1=rownames(dist_mat), dist_mat)
  dist_edge_list <- data.table::melt(dist_mat_dt, id=1, na.rm=T)
  colnames(dist_edge_list) <- c('Var1', 'Var2', 'value')
  dist_edge_list$val_rsc <- scales::rescale(
    dist_edge_list$value, to=c(0, repel_r))
  dist_edge_list$Var1 %<>% as.character
  dist_edge_list$Var2 %<>% as.character
  dist_edge_list$Var1_type <- 'sc'
  dist_edge_list$Var2_type <- 'non-sc'

  cat('Pruning graph...\n')
  dist_edge_list_sub <- dplyr::inner_join(
    dist_edge_list %>% 
      group_by(Var1) %>% 
      top_n(n=top_spot, wt=-value), 
    dist_edge_list %>% group_by(Var2) %>% 
      top_n(n=spot_n, wt=-value)) %>% 
    data.frame

  cat('Spatial Charting SC data...\n')
  sc_coord <- sc_coord_raw <- data.frame(
    id_raw=dist_edge_list_sub$Var1, id_new=make.names(dist_edge_list_sub$Var1, 
                                                      unique = T))
  sc_coord$coord_x <- 
    sc_coord_raw$coord_x <- 
    coord_df$coord_x[match(dist_edge_list_sub$Var2, rownames(coord_df))]
  sc_coord$coord_y <- 
    sc_coord_raw$coord_y <- 
    coord_df$coord_y[match(dist_edge_list_sub$Var2, rownames(coord_df))]
  ## Add noise ##
  theta <- runif(nrow(dist_edge_list_sub), 0, 2*pi)
  alpha <- sqrt(runif(nrow(dist_edge_list_sub), 0, 1))
  sc_coord$coord_x <- 
    sc_coord$coord_x + dist_edge_list_sub$val_rsc*sin(theta)*alpha
  sc_coord$coord_y <- 
    sc_coord$coord_y + dist_edge_list_sub$val_rsc*cos(theta)*alpha
  ## Point repelling ##
  cat('Repelling points...\n')
  sc_repel_input <- data.frame(
    sc_coord[, c('coord_x', 'coord_y')], repel_r=repel_r)
  sc_repel <- packcircles::circleRepelLayout(
    sc_repel_input, sizetype='radius', maxiter=repel_iter)
  sc_coord$coord_x <- sc_repel$layout$x
  sc_coord$coord_y <- sc_repel$layout$y
  return(list(sc_coord_raw, sc_coord))
}

Running the code constructs the functions of CellTrek in your environment. Now, the example of CellTrek run. Probably, other changes are necessary for specific cases. Other functions like celltrek_vis() run successfully using celltrek original library, it's mean CellTrek::celltrek_vis().

panxueling commented 3 months ago

Matrix_1.6-3 SeuratObject_4.1.3
Seurat_4.3.0 Software version problem. It was resolved after replacing the version.