Unable to subset chromatinassay by features like in vignette "Analyzing PBMC scATAC-seq" when other assays are present

[jzw1999]

In the Signac vignette "Analyzing PBMC scATAC-seq", there is a step in which we filter out those features falling in chromosome scaffolds or other sequences instead of the standard chromosomes:

peaks.keep <- seqnames(granges(pbmc)) %in% standardChromosomes(granges(pbmc))
pbmc <- pbmc[as.vector(peaks.keep), ]

I have no problem running that step when the seurat object only include a single assat "peaks". However, if the seurat object contains other assays like "RNA" (for example, the one generated in another vignette "Joint RNA and ATAC analysis: 10x multiomic"), this step will fail, and produce an error:

> pbmc2
An object of class Seurat 
204947 features across 11211 samples within 3 assays 
Active assay: ATAC (143887 features, 143887 variable features)
 2 layers present: counts, data
 2 other assays present: RNA, SCT
 2 dimensional reductions calculated: pca, lsi
> pbmc2[as.vector(peaks.keep), ]
Error: None of the features provided found in this assay

As have been mentioned in issue #315, this only happens to the features, but not the cells.

> pbmc2[, 1]
An object of class Seurat 
204947 features across 1 samples within 3 assays 
Active assay: ATAC (143887 features, 143887 variable features)
 2 layers present: counts, data
 2 other assays present: RNA, SCT
 2 dimensional reductions calculated: pca, lsi

The pbmc2 is generated as detailed in the vignette "Joint RNA and ATAC analysis: 10x multiomic", and the only modification is that I used another .h5 and fragment files from 10X genomics. I guess maybe assay RNA or SCT is causing this issue? I would really appreciate it if anyone knows the reason behind it.

I'm using the latest release of Signac (1.14.0, which should provide support for seurat v5) and Seurat 5.0.1. I would really appreciate if anyone knows the reason behind this issue.

Here is the session info:

> sessionInfo()
R version 4.2.0 (2022-04-22)
Platform: x86_64-pc-linux-gnu (64-bit)
Running under: CentOS Linux 8

Matrix products: default
BLAS/LAPACK: /software/openblas-0.3.13-el8-x86_64/lib/libopenblas_skylakexp-r0.3.13.so

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       

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

other attached packages:
 [1] SeuratDisk_0.0.0.9021             ggplot2_3.4.4                     biovizBase_1.46.0                
 [4] BSgenome.Hsapiens.UCSC.hg38_1.4.5 BSgenome_1.66.3                   rtracklayer_1.58.0               
 [7] Biostrings_2.66.0                 XVector_0.38.0                    EnsDb.Hsapiens.v86_2.99.0        
[10] ensembldb_2.22.0                  AnnotationFilter_1.22.0           GenomicFeatures_1.50.4           
[13] AnnotationDbi_1.60.2              Biobase_2.58.0                    GenomicRanges_1.50.2             
[16] GenomeInfoDb_1.34.9               IRanges_2.32.0                    S4Vectors_0.36.2                 
[19] BiocGenerics_0.44.0               Seurat_5.0.1                      SeuratObject_5.0.2               
[22] sp_1.6-0                          Signac_1.14.9001                 

loaded via a namespace (and not attached):
  [1] utf8_1.2.2                  spatstat.explore_3.2-7      reticulate_1.24             tidyselect_1.2.0           
  [5] RSQLite_2.3.1               htmlwidgets_1.6.2           grid_4.2.0                  BiocParallel_1.32.6        
  [9] Rtsne_0.16                  munsell_0.5.0               codetools_0.2-18            ica_1.0-2                  
 [13] future_1.25.0               miniUI_0.1.1.1              withr_2.5.0                 spatstat.random_3.2-3      
 [17] colorspace_2.0-3            progressr_0.10.0            filelock_1.0.2              knitr_1.42                 
 [21] rstudioapi_0.14             ROCR_1.0-11                 tensor_1.5                  listenv_0.8.0              
 [25] MatrixGenerics_1.10.0       GenomeInfoDbData_1.2.9      polyclip_1.10-0             bit64_4.0.5                
 [29] parallelly_1.31.1           vctrs_0.6.1                 generics_0.1.3              xfun_0.38                  
 [33] BiocFileCache_2.6.1         R6_2.5.1                    hdf5r_1.3.11                bitops_1.0-7               
 [37] spatstat.utils_3.0-5        cachem_1.0.6                DelayedArray_0.24.0         promises_1.2.0.1           
 [41] BiocIO_1.8.0                scales_1.3.0                nnet_7.3-17                 gtable_0.3.0               
 [45] globals_0.14.0              goftest_1.2-3               spam_2.8-0                  rlang_1.1.2                
 [49] RcppRoll_0.3.1              splines_4.2.0               lazyeval_0.2.2              dichromat_2.0-0.1          
 [53] checkmate_2.1.0             spatstat.geom_3.2-9         yaml_2.3.5                  reshape2_1.4.4             
 [57] abind_1.4-5                 backports_1.4.1             httpuv_1.6.5                Hmisc_5.1-0                
 [61] tools_4.2.0                 ellipsis_0.3.2              RColorBrewer_1.1-3          ggridges_0.5.3             
 [65] Rcpp_1.0.11                 plyr_1.8.7                  base64enc_0.1-3             progress_1.2.2             
 [69] zlibbioc_1.44.0             purrr_1.0.1                 RCurl_1.98-1.12             prettyunits_1.1.1          
 [73] rpart_4.1.16                deldir_1.0-6                pbapply_1.5-0               cowplot_1.1.1              
 [77] zoo_1.8-11                  SummarizedExperiment_1.28.0 ggrepel_0.9.3               cluster_2.1.3              
 [81] tinytex_0.38                magrittr_2.0.3              data.table_1.14.4           RSpectra_0.16-1            
 [85] scattermore_1.2             lmtest_0.9-40               RANN_2.6.1                  ProtGenerics_1.30.0        
 [89] fitdistrplus_1.1-8          matrixStats_1.1.0           evaluate_0.15               hms_1.1.3                  
 [93] patchwork_1.1.1             mime_0.12                   xtable_1.8-4                XML_3.99-0.9               
 [97] fastDummies_1.7.3           gridExtra_2.3               compiler_4.2.0              biomaRt_2.54.1             
[101] tibble_3.2.1                KernSmooth_2.23-20          crayon_1.5.1                htmltools_0.5.8.1          
[105] later_1.3.0                 Formula_1.2-4               tidyr_1.3.0                 DBI_1.1.2                  
[109] dbplyr_2.3.2                MASS_7.3-58.3               rappdirs_0.3.3              Matrix_1.6-4               
[113] cli_3.6.2                   parallel_4.2.0              dotCall64_1.0-1             igraph_2.0.3               
[117] pkgconfig_2.0.3             GenomicAlignments_1.34.1    foreign_0.8-82              plotly_4.10.0              
[121] spatstat.sparse_3.0-3       xml2_1.3.3                  VariantAnnotation_1.44.1    stringr_1.5.0              
[125] digest_0.6.29               sctransform_0.4.1           RcppAnnoy_0.0.19            spatstat.data_3.0-4        
[129] rmarkdown_2.21              leiden_0.3.10               fastmatch_1.1-3             htmlTable_2.4.0            
[133] uwot_0.1.14                 restfulr_0.0.15             curl_4.3.2                  shiny_1.7.1                
[137] Rsamtools_2.14.0            rjson_0.2.21                lifecycle_1.0.4             nlme_3.1-157               
[141] jsonlite_1.8.7              viridisLite_0.4.0           fansi_1.0.3                 pillar_1.9.0               
[145] lattice_0.20-45             KEGGREST_1.38.0             fastmap_1.2.0               httr_1.4.5                 
[149] survival_3.3-1              glue_1.6.2                  png_0.1-7                   bit_4.0.4                  
[153] stringi_1.7.6               blob_1.2.3                  RcppHNSW_0.4.1              memoise_2.0.1              
[157] dplyr_1.1.2                 irlba_2.3.3                 future.apply_1.9.0     

[jzw1999]

Actually I found a walkaround, by directly modifying the counts of assay "ATAC":

> pbmc
An object of class Seurat 
204947 features across 11211 samples within 3 assays 
Active assay: ATAC (143887 features, 143887 variable features)
 2 layers present: counts, data
 2 other assays present: RNA, SCT
 2 dimensional reductions calculated: pca, lsi
> pbmc@assays$ATAC <- subset(pbmc@assays$ATAC, features = rownames(pbmc@assays$ATAC)[as.vector(peaks.keep)])
> pbmc
An object of class Seurat 
204896 features across 11211 samples within 3 assays 
Active assay: ATAC (143836 features, 143836 variable features)
 2 layers present: counts, data
 2 other assays present: RNA, SCT
 2 dimensional reductions calculated: pca, lsi

However, I'm still interested in whether there exists another way of fixing this issue.