GetLinkedPeaks returns null even when linked peaks are present in coverage plot

[Dragonmasterx87]

Hi tim! Sorry for the bother, but when I run GetLinkedPeaks I continue to get null even though I know that the coverage plot shows links for a particular gene, why is this happening?

# insert reproducible example here
peaks_linked <- GetLinkedPeaks(master_atac, features = c("INS"), assay = "macs2", min.abs.score = 0)
> peaks_linked
NULL
> master_atac
An object of class Seurat 
737215 features across 52613 samples within 8 assays 
Active assay: macs2 (404697 features, 403130 variable features)
 7 other assays present: ATAC, RNA, prediction.score.id, predicted, chromvar, RNA_macs2, chromvar_macs2
 5 dimensional reductions calculated: lsi, umap, harmony, cca, ref.umap

> sessionInfo()
R version 4.3.2 (2023-10-31 ucrt)
Platform: x86_64-w64-mingw32/x64 (64-bit)
Running under: Windows 11 x64 (build 22621)

Matrix products: default


locale:
[1] LC_COLLATE=English_United States.utf8  LC_CTYPE=English_United States.utf8    LC_MONETARY=English_United States.utf8
[4] LC_NUMERIC=C                           LC_TIME=English_United States.utf8    

time zone: America/Chicago
tzcode source: internal

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

other attached packages:
 [1] Pando_1.0.5                       DESeq2_1.42.0                     BiocParallel_1.36.0              
 [4] cicero_1.3.9                      Gviz_1.46.1                       chromVAR_1.5.0                   
 [7] motifmatchr_1.24.0                BSgenome.Hsapiens.UCSC.hg38_1.4.5 BSgenome_1.70.1                  
[10] rtracklayer_1.62.0                BiocIO_1.12.0                     Biostrings_2.70.1                
[13] XVector_0.42.0                    TFBSTools_1.40.0                  JASPAR2020_0.99.10               
[16] qs_0.25.7                         R.utils_2.12.3                    R.oo_1.25.0                      
[19] R.methodsS3_1.8.2                 devtools_2.4.5                    usethis_2.2.2                    
[22] ggVennDiagram_1.5.0               ggvenn_0.1.10                     DropletUtils_1.22.0              
[25] Nebulosa_1.12.0                   circlize_0.4.15                   ComplexHeatmap_2.18.0            
[28] viridis_0.6.4                     viridisLite_0.4.2                 EnrichmentBrowser_2.32.0         
[31] graph_1.80.0                      escape_1.12.0                     dittoSeq_1.14.0                  
[34] DOSE_3.28.2                       clusterProfiler_4.10.0            MeSHDbi_1.38.0                   
[37] AnnotationHub_3.10.0              BiocFileCache_2.10.1              dbplyr_2.4.0                     
[40] org.Hs.eg.db_3.18.0               GOSemSim_2.28.1                   glmGamPoi_1.14.0                 
[43] EnhancedVolcano_1.20.0            DoubletFinder_2.0.4               future_1.33.1                    
[46] patchwork_1.2.0                   clustree_0.5.1                    ggraph_2.1.0                     
[49] plotly_4.10.4                     EnsDb.Hsapiens.v86_2.99.0         ensembldb_2.26.0                 
[52] AnnotationFilter_1.26.0           GenomicFeatures_1.54.1            AnnotationDbi_1.64.1             
[55] scDblFinder_1.16.0                Signac_1.11.0                     harmony_1.2.0                    
[58] monocle3_1.3.4                    SingleCellExperiment_1.24.0       SummarizedExperiment_1.32.0      
[61] GenomicRanges_1.54.1              GenomeInfoDb_1.38.5               IRanges_2.36.0                   
[64] S4Vectors_0.40.2                  MatrixGenerics_1.14.0             matrixStats_1.2.0                
[67] Biobase_2.62.0                    BiocGenerics_0.48.1               SeuratObject_4.1.4               
[70] Seurat_4.3.0                      reticulate_1.34.0                 data.table_1.14.10               
[73] lubridate_1.9.3                   forcats_1.0.0                     purrr_1.0.2                      
[76] readr_2.1.5                       tidyr_1.3.0                       tibble_3.2.1                     
[79] tidyverse_2.0.0                   dplyr_1.1.4                       ggridges_0.5.5                   
[82] Matrix_1.6-5                      cowplot_1.1.2                     Rcpp_1.0.12                      
[85] SoupX_1.6.2                       hdf5r_1.3.9                       stringr_1.5.1                    
[88] leiden_0.4.3.1                    ggrepel_0.9.5                     ggplot2_3.4.4                    

loaded via a namespace (and not attached):
  [1] igraph_1.6.0                  Formula_1.2-5                 ica_1.0-3                     scater_1.30.1                
  [5] maps_3.4.2                    zlibbioc_1.48.0               tidyselect_1.2.0              bit_4.0.5                    
  [9] doParallel_1.0.17             clue_0.3-65                   lattice_0.22-5                rjson_0.2.21                 
 [13] blob_1.2.4                    urlchecker_1.0.1              S4Arrays_1.2.0                dichromat_2.0-0.1            
 [17] parallel_4.3.2                seqLogo_1.68.0                png_0.1-8                     cli_3.6.2                    
 [21] ggplotify_0.1.2               ProtGenerics_1.34.0           goftest_1.2-3                 bluster_1.12.0               
 [25] BiocNeighbors_1.20.2          uwot_0.1.16                   shadowtext_0.1.3              curl_5.2.0                   
 [29] evaluate_0.23                 mime_0.12                     tidytree_0.4.6                stringi_1.8.3                
 [33] backports_1.4.1               XML_3.99-0.16                 httpuv_1.6.13                 magrittr_2.0.3               
 [37] rappdirs_0.3.3                splines_4.3.2                 RcppRoll_0.3.0                mclust_6.0.1                 
 [41] RApiSerialize_0.1.2           jpeg_0.1-10                   DT_0.31                       sctransform_0.4.1            
 [45] ggbeeswarm_0.7.2              sessioninfo_1.2.2             DBI_1.2.1                     terra_1.7-65                 
 [49] HDF5Array_1.30.0              withr_3.0.0                   rprojroot_2.0.4               enrichplot_1.22.0            
 [53] xgboost_1.7.6.1               lmtest_0.9-40                 GSEABase_1.64.0               tidygraph_1.3.0              
 [57] BiocManager_1.30.22           htmlwidgets_1.6.4             fs_1.6.3                      biomaRt_2.58.0               
 [61] labeling_0.4.3                SparseArray_1.2.3             annotate_1.80.0               VariantAnnotation_1.48.1     
 [65] zoo_1.8-12                    knitr_1.45                    TFMPvalue_0.0.9               timechange_0.3.0             
 [69] foreach_1.5.2                 ggpointdensity_0.1.0          fansi_1.0.6                   caTools_1.18.2               
 [73] ggtree_3.10.0                 rhdf5_2.46.1                  poweRlaw_0.70.6               irlba_2.3.5.1                
 [77] ggrastr_1.0.2                 gridGraphics_0.5-1            ellipsis_0.3.2                lazyeval_0.2.2               
 [81] yaml_2.3.8                    survival_3.5-7                scattermore_1.2               BiocVersion_3.18.1           
 [85] crayon_1.5.2                  RcppAnnoy_0.0.21              mapproj_1.2.11                RColorBrewer_1.1-3           
 [89] progressr_0.14.0              tweenr_2.0.2                  later_1.3.2                   Rgraphviz_2.46.0             
 [93] base64enc_0.1-3               codetools_0.2-19              GlobalOptions_0.1.2           profvis_0.3.8                
 [97] KEGGREST_1.42.0               Rtsne_0.17                    shape_1.4.6                   limma_3.58.1                 
[101] Rsamtools_2.18.0              filelock_1.0.3                foreign_0.8-86                pkgconfig_2.0.3              
[105] KEGGgraph_1.62.0              xml2_1.3.6                    GenomicAlignments_1.38.2      aplot_0.2.2                  
[109] biovizBase_1.50.0             spatstat.sparse_3.0-3         ape_5.7-1                     xtable_1.8-4                 
[113] interp_1.1-5                  plyr_1.8.9                    httr_1.4.7                    tools_4.3.2                  
[117] globals_0.16.2                pkgbuild_1.4.3                checkmate_2.3.1               htmlTable_2.4.2              
[121] beeswarm_0.4.0                broom_1.0.5                   nlme_3.1-164                  HDO.db_0.99.1                
[125] lme4_1.1-35.1                 digest_0.6.34                 farver_2.1.1                  tzdb_0.4.0                   
[129] reshape2_1.4.4                ks_1.14.2                     yulab.utils_0.1.3             rpart_4.1.23                 
[133] DirichletMultinomial_1.44.0   glue_1.7.0                    cachem_1.0.8                  polyclip_1.10-6              
[137] Hmisc_5.1-1                   generics_0.1.3                mvtnorm_1.2-4                 parallelly_1.36.0            
[141] pkgload_1.3.4                 statmod_1.5.0                 here_1.0.1                    ScaledMatrix_1.10.0          
[145] minqa_1.2.6                   pbapply_1.7-2                 nabor_0.5.0                   spam_2.10-0                  
[149] gson_0.1.0                    dqrng_0.3.2                   utf8_1.2.4                    gtools_3.9.5                 
[153] graphlayouts_1.1.0            gridExtra_2.3                 shiny_1.8.0                   GSVA_1.50.0                  
[157] GenomeInfoDbData_1.2.11       pals_1.8                      rhdf5filters_1.14.1           RCurl_1.98-1.14              
[161] memoise_2.0.1                 rmarkdown_2.25                pheatmap_1.0.12               scales_1.3.0                 
[165] RANN_2.6.1                    stringfish_0.16.0             spatstat.data_3.0-4           rstudioapi_0.15.0            
[169] cluster_2.1.6                 msigdbr_7.5.1                 spatstat.utils_3.0-4          hms_1.1.3                    
[173] fitdistrplus_1.1-11           munsell_0.5.0                 colorspace_2.1-0              rlang_1.1.3                  
[177] DelayedMatrixStats_1.24.0     sparseMatrixStats_1.14.0      dotCall64_1.1-1               ggforce_0.4.1                
[181] scuttle_1.12.0                xfun_0.41                     CNEr_1.38.0                   remotes_2.4.2.1              
[185] iterators_1.0.14              abind_1.4-5                   interactiveDisplayBase_1.40.0 treeio_1.26.0                
[189] Rhdf5lib_1.24.1               bitops_1.0-7                  promises_1.2.1                scatterpie_0.2.1             
[193] RSQLite_2.3.4                 qvalue_2.34.0                 fgsea_1.28.0                  DelayedArray_0.28.0          
[197] GO.db_3.18.0                  compiler_4.3.2                prettyunits_1.2.0             boot_1.3-28.1                
[201] beachmat_2.18.0               listenv_0.9.0                 edgeR_4.0.10                  BiocSingular_1.18.0          
[205] tensor_1.5                    MASS_7.3-60.0.1               progress_1.2.3                UCell_2.6.2                  
[209] babelgene_22.9                spatstat.random_3.2-2         R6_2.5.1                      fastmap_1.1.1                
[213] fastmatch_1.1-4               vipor_0.4.7                   ROCR_1.0-11                   nnet_7.3-19                  
[217] rsvd_1.0.5                    gtable_0.3.4                  KernSmooth_2.23-22            latticeExtra_0.6-30          
[221] miniUI_0.1.1.1                deldir_2.0-2                  htmltools_0.5.7               RcppParallel_5.1.7           
[225] bit64_4.0.5                   spatstat.explore_3.2-5        lifecycle_1.0.4               nloptr_2.0.3                 
[229] restfulr_0.0.15               vctrs_0.6.5                   VGAM_1.1-9                    spatstat.geom_3.2-7          
[233] scran_1.30.1                  ggfun_0.1.4                   sp_2.1-2                      pracma_2.4.4                 
[237] future.apply_1.11.1           pillar_1.9.0                  metapod_1.10.1                locfit_1.5-9.8               
[241] jsonlite_1.8.8                GetoptLong_1.0.5 

My goal is to identify Motifs within peaks linked to a gene, say INS, and compare across samples the accessibility of specific motifs in proximity to genes say INS. Links are present for individually calculated cellranger and macs2 outputs:

> master_atac@assays[["macs2"]]@links
GRanges object with 20159992 ranges and 2 metadata columns:
             seqnames              ranges strand |       score     group
                <Rle>           <IRanges>  <Rle> |   <numeric> <numeric>
         [1]     chr1  99784242-100037921      * | 0.000632882      1295
         [2]     chr1  99836316-100037921      * | 0.046875528      1295
         [3]     chr1  99850045-100037921      * | 0.099272697      1289
         [4]     chr1  99969902-100037921      * | 0.111595901      1312
         [5]     chr1 100037921-100065775      * | 0.016352345      1308
         ...      ...                 ...    ... .         ...       ...
  [20159988]     chrY     7465401-7718738      * |   0.0028103        NA
  [20159989]     chrY     8712378-8712942      * |   0.4171982        NA
  [20159990]     chrY     8712378-8712942      * |   0.4171982        NA
  [20159991]     chrY     8780958-8781461      * |   0.2446885        NA
  [20159992]     chrY     8780958-8781461      * |   0.2446885        NA
  -------
  seqinfo: 24 sequences from an unspecified genome; no seqlengths

Cheers for the help!

🐉

[timoast]

Hi @Dragonmasterx87!

Are these links from cicero? If so, cicero clusters groups of peaks together, rather than linking peaks to genes. The GetLinkedPeaks() function is meant for extracting information generated by the LinkPeaks() function that finds the correlation between peak accessibility and gene expression (although I now realize this is not well documented, I will update it soon)

[timoast]

Documentation is now updated