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

Activity

1. tsoruu commented on Jan 9, 2024

   

   tsoruu

   on Jan 9, 2024

   Author

   More actions

   Got some progress on this issue, I've tried converting the object to a similar format to v3 using satijalab/seurat#7409 (comment) 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?
   
2. Mapleljs commented on Jan 23, 2024

   

   Mapleljs

   on Jan 23, 2024

   More actions

   I am facing the same question as you. Have you known how to solve it?
3. xieaq commented on Mar 1, 2024

   

   xieaq

   on Mar 1, 2024

   More actions

   Same question :(
4. AndoneguiElguera commented on Apr 3, 2024

   

   AndoneguiElguera

   on Apr 3, 2024

   Last edited by AndoneguiElguera

   More actions

   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().

   👍React with 👍3Reacted by YiDing, diego4-2 and Zerun Song
5. panxueling commented on Jun 12, 2024

   

   panxueling

   on Jun 12, 2024

   Last edited by panxueling

   More actions

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

   👍React with 👍2Reacted by Yann Maugé and sixi
6. TangentDeng commented on Nov 12, 2024

   

   TangentDeng

   on Nov 12, 2024

   More actions

   @panxueling what version are you used?
7. NIKO6666666 commented on Nov 28, 2024

   

   NIKO6666666

   on Nov 28, 2024

   More actions

   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().

   thanks, my dear friends.

   👍React with 👍1Reacted by Zerun Song
8. YannMauge commented on Dec 10, 2024

   

   YannMauge

   on Dec 10, 2024

   More actions

   @panxueling what version are you used?
   I have been running fine with SeuratObject_4.1.4 / Seurat_4.4.0
   SeuratObject_4.1.3 / Seurat_4.3.0 should also work as mentionned by @panxueling if you don't want to change the functions for Seurat 5.