Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

PrepSCTIntegration Error in scale.data[anchor.features, ] : subscript out of bounds #8216

Closed
LukeJ89 opened this issue Dec 21, 2023 · 5 comments

Comments

@LukeJ89
Copy link

LukeJ89 commented Dec 21, 2023

Hi,

I'm fairly new to R and Seurat, so I hope this will be an easy fix. I've a list of 7 Seurat objects from 10X spatial datasets, and ran SCTransform on them. The list is "Prep_Healthy". Running this code:

Prep_Healthy<- SCTransform(Prep_Healthy, vst.flavor = "v2", assay="Spatial")
features <- SelectIntegrationFeatures(object.list = Prep_healthy, nfeatures = 3000)
Prep_healthy <- PrepSCTIntegration(object.list = Prep_healthy, anchor.features = features)

But running this I get the error code "| 0 % ~calculating Error in scale.data[anchor.features, ] : subscript out of bounds"
If I set the nfeatures lower, then it works, the highest is nfeatures = 690 after which the error occurs. So it seems that the list for the anchor.features is too long, not sure if it's a memory thing but I've tried it on 2 computers.

I also tried merging and splitting the list as in #3085 and also tried to use intersect in case some names were missing across the datasets as in #6707. I'm at a bit of a loss, any help would be greatly appreciated.

"
Healthy_Samples <- list(HEALTHY.Male.s1, HEALTHY.Male.s2,
HEALTHY.Female1.s1, HEALTHY.Female1.s2.clean,
HEALTHY.Female1.s3, HEALTHY.Female2.s1, HEALTHY.Female2.s2)

#Filter for counts greater than 200
i <- 1
while(i <= length(Healthy_Samples)){
filtered_data <- st_filter_by_genes(st.data = Healthy_Samples[[i]],x = 200)
Healthy_Samples[[i]] <- filtered_data
i <- i+1
}

st_SCT <- function(data.list){
for (i in 1:length(data.list)) {
data.list[[i]] <- SCTransform(data.list[[i]], vst.flavor = "v2", verbose = T,assay="Spatial")}
data.list[[i]] <- UpdateSeuratObject(data.list[[i]])
return(data.list)
}

Prep_healthy <- st_SCT(Healthy_Samples)

features <- SelectIntegrationFeatures(object.list = Prep_healthy, nfeatures = 3000)
Prep_healthy <- PrepSCTIntegration(object.list = Prep_healthy, anchor.features = features)

@saketkc
Copy link
Collaborator

saketkc commented Dec 21, 2023

features <- SelectIntegrationFeatures(object.list = Prep_healthy, nfeatures = 3000)

expects that object.list to be list. Assuming this was an issue while pasting, can you paste your entire code here along with sessionInfo()

@LukeJ89
Copy link
Author

LukeJ89 commented Dec 21, 2023

st_combine <- function(data.list,ndim,res,seed=100){
  set.seed(seed)
  for (i in 1:length(data.list)) {
    data.list[[i]] <- SCTransform(data.list[[i]], vst.flavor = "v2", verbose = T,assay="Spatial")
  }
  features <- SelectIntegrationFeatures(object.list =  data.list, nfeatures = 3000)
  data.list <- PrepSCTIntegration(object.list = data.list, anchor.features = features)
  data.anchors <- FindIntegrationAnchors(object.list = data.list, normalization.method = "SCT", 
                                         anchor.features = features)
  data.combined.sct <- IntegrateData(anchorset = data.anchors, normalization.method = "SCT")
  data.combined.sct <- RunPCA(data.combined.sct, verbose = FALSE)
  data.combined.sct <- FindNeighbors(data.combined.sct, dims = 1:ndim)
  data.combined.sct <- FindClusters(data.combined.sct, verbose = FALSE,resolution = res)
  data.combined.sct <- RunUMAP(data.combined.sct, dims = 1:ndim)
  return(data.combined.sct)
}


Prep_healthy <- c(HEALTHY.Male.s1, HEALTHY.Male.s2,
                        HEALTHY.Female1.s1, HEALTHY.Female1.s2.clean,
                        HEALTHY.Female1.s3, HEALTHY.Female2.s1, HEALTHY.Female2.s2)

#Filter for counts greater than 200
i <- 1
while(i <= length(Prep_healthy)){
  filtered_data <- st_filter_by_genes(st.data = Prep_healthy[[i]],x = 200)
  Prep_healthy[[i]] <- filtered_data
  i <- i+1
}

Prep_healthy <- st_combine(Prep_healthy, ndim = 20, res = 0.6)

Each name in the "Prep_healthy" list is a Seurat object

  |                                                  | 0 % ~calculating  Error in scale.data[anchor.features, ] : subscript out of bounds
In addition: There were 50 or more warnings (use warnings() to see the first 50)


> sessionInfo()
R version 4.2.3 (2023-03-15 ucrt)
Platform: x86_64-w64-mingw32/x64 (64-bit)
Running under: Windows 10 x64 (build 22635)

Matrix products: default

locale:
[1] LC_COLLATE=English_United Kingdom.utf8  LC_CTYPE=English_United Kingdom.utf8   
[3] LC_MONETARY=English_United Kingdom.utf8 LC_NUMERIC=C                           
[5] LC_TIME=English_United Kingdom.utf8    

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

other attached packages:
 [1] harmony_1.1.0      Rcpp_1.0.11        ggsci_3.0.0        scales_1.3.0       pheatmap_1.0.12   
 [6] RColorBrewer_1.1-3 cowplot_1.1.2      Seurat_5.0.1.9001  SeuratObject_5.0.1 sp_2.1-1          
[11] lubridate_1.9.3    forcats_1.0.0      stringr_1.5.1      dplyr_1.1.4        purrr_1.0.2       
[16] readr_2.1.4        tidyr_1.3.0        tibble_3.2.1       ggplot2_3.4.4      tidyverse_2.0.0   

loaded via a namespace (and not attached):
  [1] spam_2.10-0                 plyr_1.8.9                  igraph_1.5.1               
  [4] lazyeval_0.2.2              splines_4.2.3               RcppHNSW_0.5.0             
  [7] listenv_0.9.0               scattermore_1.2             GenomeInfoDb_1.34.9        
 [10] usethis_2.2.2               digest_0.6.33               htmltools_0.5.7            
 [13] fansi_1.0.5                 magrittr_2.0.3              memoise_2.0.1              
 [16] tensor_1.5                  cluster_2.1.4               ROCR_1.0-11                
 [19] tzdb_0.4.0                  remotes_2.4.2.1             globals_0.16.2             
 [22] matrixStats_1.1.0           timechange_0.2.0            spatstat.sparse_3.0-3      
 [25] prettyunits_1.2.0           colorspace_2.1-0            ggrepel_0.9.4              
 [28] RCurl_1.98-1.13             callr_3.7.3                 crayon_1.5.2               
 [31] jsonlite_1.8.7              progressr_0.14.0            spatstat.data_3.0-3        
 [34] survival_3.5-7              zoo_1.8-12                  glue_1.6.2                 
 [37] polyclip_1.10-6             gtable_0.3.4                zlibbioc_1.44.0            
 [40] XVector_0.38.0              leiden_0.4.3.1              DelayedArray_0.24.0        
 [43] pkgbuild_1.4.2              future.apply_1.11.0         BiocGenerics_0.44.0        
 [46] abind_1.4-5                 spatstat.random_3.2-1       miniUI_0.1.1.1             
 [49] viridisLite_0.4.2           xtable_1.8-4                reticulate_1.34.0          
 [52] dotCall64_1.1-0             stats4_4.2.3                profvis_0.3.8              
 [55] htmlwidgets_1.6.4           httr_1.4.7                  ellipsis_0.3.2             
 [58] ica_1.0-3                   urlchecker_1.0.1            pkgconfig_2.0.3            
 [61] farver_2.1.1                uwot_0.1.16                 deldir_2.0-2               
 [64] utf8_1.2.4                  tidyselect_1.2.0            labeling_0.4.3             
 [67] rlang_1.1.2                 reshape2_1.4.4              later_1.3.1                
 [70] munsell_0.5.0               tools_4.2.3                 cachem_1.0.8               
 [73] cli_3.6.1                   generics_0.1.3              devtools_2.4.5             
 [76] ggridges_0.5.5              fastmap_1.1.1               goftest_1.2-3              
 [79] processx_3.8.2              fs_1.6.3                    fitdistrplus_1.1-11        
 [82] RANN_2.6.1                  sparseMatrixStats_1.10.0    pbapply_1.7-2              
 [85] future_1.33.0               nlme_3.1-163                mime_0.12                  
 [88] compiler_4.2.3              rstudioapi_0.15.0           plotly_4.10.3              
 [91] png_0.1-8                   spatstat.utils_3.0-4        glmGamPoi_1.10.2           
 [94] stringi_1.8.2               ps_1.7.5                    RSpectra_0.16-1            
 [97] lattice_0.22-5              Matrix_1.6-3                vctrs_0.6.4                
[100] pillar_1.9.0                lifecycle_1.0.4             BiocManager_1.30.22        
[103] spatstat.geom_3.2-7         lmtest_0.9-40               RcppAnnoy_0.0.21           
[106] bitops_1.0-7                data.table_1.14.8           irlba_2.3.5.1              
[109] GenomicRanges_1.50.2        httpuv_1.6.12               patchwork_1.1.3            
[112] R6_2.5.1                    promises_1.2.1              KernSmooth_2.23-22         
[115] gridExtra_2.3               IRanges_2.32.0              parallelly_1.36.0          
[118] sessioninfo_1.2.2           codetools_0.2-19            fastDummies_1.7.3          
[121] MASS_7.3-60                 pkgload_1.3.3               SummarizedExperiment_1.28.0
[124] withr_2.5.2                 sctransform_0.4.1           GenomeInfoDbData_1.2.9     
[127] S4Vectors_0.36.2            parallel_4.2.3              hms_1.1.3                  
[130] grid_4.2.3                  DelayedMatrixStats_1.20.0   MatrixGenerics_1.10.0      
[133] Rtsne_0.16                  spatstat.explore_3.2-5      Biobase_2.58.0             
[136] shiny_1.8.0   

@saketkc
Copy link
Collaborator

saketkc commented Dec 21, 2023

Ahh, I see. This is related to #8203. My suggestion for now would be to create a "RNA" assay, by copying the "Spatial" assay and doing everything on the "RNA" assay.

object[["RNA"]] <- object[["Spatial"]]

@LukeJ89
Copy link
Author

LukeJ89 commented Dec 21, 2023

Perfect! That works, thank you so much!

@LukeJ89 LukeJ89 closed this as completed Dec 21, 2023
@FerrenaAlexander
Copy link

Hi all, thanks for all of your amazing hard work!
Just wanted to report something maybe related to this, I was able to get around the following error in the GetResidual function: "Error in .subscript.2ary(x, , j, drop = drop) : subscript out of bounds".

By setting object[["RNA"]] <- object[["Spatial"]] , at two specific times:

  1. before running SCTransform on each indiviudal object
  2. before running PrepSCTIntegration

I was using GetResidual after integration, to plot genes that were not in the scale.data layer after integration. Not sure if the second step was really necessary. Doing the first of these fixed the GetResidual error in integration, but the second fixed only the PrepSCTIntegration error.

Also not sure whether this may be related to #8561.

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

No branches or pull requests

3 participants