Pando
Pando copied to clipboard
Published
20 hours ago •
quadbio
quadbio
Error in infer_grn()
Dear Pndo developers,
I was wondering if you could help me to debug my code - I am just following your tutorial and applying it on my data below, however I ran into an error irrespective of which setup I use:
pbmc10
An object of class Seurat
601236 features across 116471 samples within 4 assays
Active assay: peaks (349688 features, 0 variable features)
3 other assays present: RNA, ATAC, SCT
7 dimensional reductions calculated: pca, harmony.pca, umap.rnah, lsi, harmony.lsi, umap.atach, wnn.umap
#---initialize GRN object
pbmc10 <- initiate_grn(pbmc10)
library(BSgenome.Hsapiens.UCSC.hg38)
data(motifs)
pbmc10 <- find_motifs(
pbmc10,
pfm = motifs,
genome = BSgenome.Hsapiens.UCSC.hg38
)
pbmc10 <- infer_grn(
pbmc10#,
#peak_to_gene_method = 'GREAT',
#method = 'brms'
)
Selecting candidate regulatory regions near genes
Error in .subset2(x, i, exact = exact): attempt to select less than one element in integerOneIndex
Traceback:
1. infer_grn(pbmc10)
2. infer_grn.SeuratPlus(pbmc10)
3. fit_grn_models(object = object, genes = genes, network_name = network_name,
. peak_to_gene_method = peak_to_gene_method, upstream = upstream,
. downstream = downstream, extend = extend, only_tss = only_tss,
. parallel = parallel, tf_cor = tf_cor, peak_cor = peak_cor,
. aggregate_rna_col = aggregate_rna_col, aggregate_peaks_col = aggregate_peaks_col,
. method = method, alpha = alpha, family = family, interaction_term = interaction_term,
. adjust_method = adjust_method, scale = scale, verbose = verbose,
. ...)
4. fit_grn_models.SeuratPlus(object = object, genes = genes, network_name = network_name,
. peak_to_gene_method = peak_to_gene_method, upstream = upstream,
. downstream = downstream, extend = extend, only_tss = only_tss,
. parallel = parallel, tf_cor = tf_cor, peak_cor = peak_cor,
. aggregate_rna_col = aggregate_rna_col, aggregate_peaks_col = aggregate_peaks_col,
. method = method, alpha = alpha, family = family, interaction_term = interaction_term,
. adjust_method = adjust_method, scale = scale, verbose = verbose,
. ...)
5. aggregate_matrix(t(peaks_near_gene), groups = colnames(peaks_near_gene),
. fun = "sum")
6. fast_aggregate(x = x, groupings = groups, fun = fun)
7. data.frame(interaction(groupings2, sep = "_"))
8. interaction(groupings2, sep = "_")
9. as.factor(args[[i]])
10. is.factor(x)
11. args[[i]]
12. `[[.data.frame`(args, i)
13. (function(x, i, exact) if (is.matrix(i)) as.matrix(x)[[i]] else .subset2(x,
. i, exact = exact))(x, ..., exact = exact)
I am not sure what could be the problem, I tried to DietSeurat, so it only contains RNA and peaks or RNA and ATAC assay but I still got the same problem.
Thank you so much for your help.
Best regards, Tereza
>> sessionInfo()
R version 4.1.1 (2021-08-10)
Platform: x86_64-conda-linux-gnu (64-bit)
Running under: Ubuntu 21.04
Matrix products: default
locale:
[1] C
attached base packages:
[1] stats4 stats graphics grDevices utils datasets methods
[8] base
other attached packages:
[1] BSgenome.Hsapiens.UCSC.hg38_1.4.4 BSgenome_1.62.0
[3] rtracklayer_1.54.0 Biostrings_2.62.0
[5] XVector_0.34.0 Pando_1.0.3
[7] regioneR_1.26.1 GenomicRanges_1.46.1
[9] GenomeInfoDb_1.30.1 IRanges_2.28.0
[11] S4Vectors_0.32.4 BiocGenerics_0.40.0
[13] dplyr_1.0.10 Signac_1.9.0
[15] SeuratObject_4.1.3 Seurat_4.3.0
[17] dreamlet_0.0.39 variancePartition_1.25.13
[19] BiocParallel_1.28.3 limma_3.50.3
[21] ggplot2_3.4.0
loaded via a namespace (and not attached):
[1] pbdZMQ_0.3-8 scattermore_0.8
[3] mashr_0.2.69 R.methodsS3_1.8.2
[5] tidyr_1.2.1 bit64_4.0.5
[7] R.utils_2.12.2 irlba_2.3.5.1
[9] DelayedArray_0.20.0 data.table_1.14.6
[11] TFBSTools_1.32.0 KEGGREST_1.34.0
[13] RCurl_1.98-1.9 doParallel_1.0.17
[15] generics_0.1.3 RhpcBLASctl_0.21-247.1
[17] cowplot_1.1.1 RSQLite_2.2.20
[19] RANN_2.6.1 future_1.30.0
[21] ggpointdensity_0.1.0 tzdb_0.3.0
[23] bit_4.0.5 spatstat.data_3.0-0
[25] httpuv_1.6.7 SummarizedExperiment_1.24.0
[27] assertthat_0.2.1 DirichletMultinomial_1.36.0
[29] viridis_0.6.2 hms_1.1.2
[31] babelgene_22.9 evaluate_0.19
[33] promises_1.2.0.1 fansi_1.0.3
[35] restfulr_0.0.15 progress_1.2.2
[37] caTools_1.18.2 Rgraphviz_2.38.0
[39] igraph_1.3.5 DBI_1.1.3
[41] htmlwidgets_1.6.0 spatstat.geom_3.0-3
[43] purrr_1.0.0 ellipsis_0.3.2
[45] backports_1.4.1 annotate_1.72.0
[47] aod_1.3.2 deldir_1.0-6
[49] sparseMatrixStats_1.6.0 MatrixGenerics_1.6.0
[51] vctrs_0.5.1 SingleCellExperiment_1.16.0
[53] Biobase_2.54.0 ROCR_1.0-11
[55] abind_1.4-5 cachem_1.0.6
[57] withr_2.5.0 ggforce_0.4.1
[59] progressr_0.12.0 sctransform_0.3.5
[61] GenomicAlignments_1.30.0 prettyunits_1.1.1
[63] goftest_1.2-3 cluster_2.1.4
[65] seqLogo_1.60.0 IRdisplay_1.1
[67] lazyeval_0.2.2 crayon_1.5.2
[69] spatstat.explore_3.0-5 edgeR_3.36.0
[71] pkgconfig_2.0.3 tweenr_2.0.2
[73] nlme_3.1-161 pals_1.7
[75] rlang_1.0.6 globals_0.16.2
[77] lifecycle_1.0.3 miniUI_0.1.1.1
[79] dichromat_2.0-0.1 invgamma_1.1
[81] polyclip_1.10-4 matrixStats_0.63.0
[83] lmtest_0.9-40 graph_1.72.0
[85] Matrix_1.5-3 IRkernel_1.3.1
[87] ashr_2.2-54 Rhdf5lib_1.16.0
[89] boot_1.3-28.1 zoo_1.8-11
[91] base64enc_0.1-3 ggridges_0.5.4
[93] GlobalOptions_0.1.2 png_0.1-8
[95] viridisLite_0.4.1 rjson_0.2.21
[97] bitops_1.0-7 R.oo_1.25.0
[99] KernSmooth_2.23-20 rhdf5filters_1.6.0
[101] EnrichmentBrowser_2.25.3 blob_1.2.3
[103] DelayedMatrixStats_1.16.0 shape_1.4.6
[105] mixsqp_0.3-48 stringr_1.5.0
[107] SQUAREM_2021.1 parallelly_1.33.0
[109] spatstat.random_3.0-1 readr_2.1.3
[111] CNEr_1.30.0 rmeta_3.0
[113] scales_1.2.1 memoise_2.0.1
[115] GSEABase_1.56.0 magrittr_2.0.3
[117] plyr_1.8.8 ica_1.0-3
[119] gplots_3.1.3 zlibbioc_1.40.0
[121] compiler_4.1.1 BiocIO_1.4.0
[123] RColorBrewer_1.1-3 clue_0.3-63
[125] KEGGgraph_1.54.0 lme4_1.1-31
[127] fitdistrplus_1.1-8 Rsamtools_2.10.0
[129] cli_3.5.0 listenv_0.9.0
[131] patchwork_1.1.2 pbapply_1.6-0
[133] MASS_7.3-58.1 tidyselect_1.2.0
[135] stringi_1.7.8 zenith_0.99.1
[137] yaml_2.3.6 locfit_1.5-9.7
[139] ggrepel_0.9.2 grid_4.1.1
[141] fastmatch_1.1-3 tools_4.1.1
[143] future.apply_1.10.0 parallel_4.1.1
[145] circlize_0.4.15 uuid_1.1-0
[147] TFMPvalue_0.0.9 foreach_1.5.2
[149] gridExtra_2.3 farver_2.1.1
[151] Rtsne_0.16 ggraph_2.1.0
[153] RcppZiggurat_0.1.6 digest_0.6.31
[155] pracma_2.4.2 shiny_1.7.4
[157] motifmatchr_1.16.0 Rcpp_1.0.9
[159] broom_1.0.2 later_1.3.0
[161] RcppAnnoy_0.0.20 httr_1.4.4
[163] AnnotationDbi_1.56.2 ComplexHeatmap_2.10.0
[165] Rdpack_2.4 colorspace_2.0-3
[167] XML_3.99-0.13 tensor_1.5
[169] reticulate_1.26 truncnorm_1.0-8
[171] splines_4.1.1 uwot_0.1.14
[173] RcppRoll_0.3.0 spatstat.utils_3.0-1
[175] graphlayouts_0.8.4 sp_1.5-1
[177] mapproj_1.2.11 plotly_4.10.1
[179] xtable_1.8-4 poweRlaw_0.70.6
[181] jsonlite_1.8.4 nloptr_2.0.3
[183] tidygraph_1.2.2 Rfast_2.0.6
[185] R6_2.5.1 pillar_1.8.1
[187] htmltools_0.5.4 mime_0.12
[189] glue_1.6.2 fastmap_1.1.0
[191] minqa_1.2.5 codetools_0.2-18
[193] maps_3.4.1 mvtnorm_1.1-3
[195] utf8_1.2.2 lattice_0.20-45
[197] spatstat.sparse_3.0-0 tibble_3.1.8
[199] pbkrtest_0.5.1 leiden_0.4.3
[201] gtools_3.9.4 GO.db_3.14.0
[203] survival_3.4-0 repr_1.1.4
[205] munsell_0.5.0 GetoptLong_1.0.5
[207] rhdf5_2.38.1 GenomeInfoDbData_1.2.7
[209] iterators_1.0.14 HDF5Array_1.22.1
[211] reshape2_1.4.4 gtable_0.3.1
[213] msigdbr_7.5.1 rbibutils_2.2.11
Hi @TerezaClarence , I'm not fully sure what causes the issue here, are you using a publically PBMC available dataset? If you point me to where you have it from then I can try to reproduce the error
Dear @joschif I am using my multiome dataset, data were loaded from cellranger 10x multiome and preprocessed with seurat and signac. I have no idea why it's giving me the error - I tried to repeatedly just diet the data so it contains RNA and ATC-assay, but it still gives me the same error. Changing the DefaultAssay also didn't help.
Best, Tereza
Hello @TerezaClarence! Just to make sure, the RNA assay does have variable features, correct?
Hello,
I used the scMultiome analysis vignette on our unpublished data in which in Step3 a gene regulatory network reconstruction is performed using pando.
When I run infer_grn I get the same error message as mentioned above:
attempt to select less than one element in integerOneIndex
13: (function(x, i, exact) if (is.matrix(i)) as.matrix(x)[[i]] else .subset2(x,
i, exact = exact))(x, ..., exact = exact)
12: `[[.data.frame`(args, i)
11: args[[i]]
10: is.factor(x)
9: as.factor(args[[i]])
8: interaction(groupings2, sep = "_")
7: data.frame(interaction(groupings2, sep = "_"))
6: fast_aggregate(x = x, groupings = groups, fun = fun)
5: aggregate_matrix(t(peaks_near_gene), groups = colnames(peaks_near_gene),
fun = "sum")
4: fit_grn_models.SeuratPlus(object = object, genes = genes, network_name = network_name,
peak_to_gene_method = peak_to_gene_method, upstream = upstream,
downstream = downstream, extend = extend, only_tss = only_tss,
parallel = parallel, tf_cor = tf_cor, peak_cor = peak_cor,
aggregate_rna_col = aggregate_rna_col, aggregate_peaks_col = aggregate_peaks_col,
method = method, alpha = alpha, family = family, interaction_term = interaction_term,
adjust_method = adjust_method, scale = scale, verbose = verbose,
...)
3: fit_grn_models(object = object, genes = genes, network_name = network_name,
peak_to_gene_method = peak_to_gene_method, upstream = upstream,
downstream = downstream, extend = extend, only_tss = only_tss,
parallel = parallel, tf_cor = tf_cor, peak_cor = peak_cor,
aggregate_rna_col = aggregate_rna_col, aggregate_peaks_col = aggregate_peaks_col,
method = method, alpha = alpha, family = family, interaction_term = interaction_term,
adjust_method = adjust_method, scale = scale, verbose = verbose,
...)
2: infer_grn.SeuratPlus(seurat, parallel = F, tf_cor = 0.05, method = "glm",
family = "gaussian", scale = F, verbose = T)
1: infer_grn(seurat, parallel = F, tf_cor = 0.05, method = "glm",
family = "gaussian", scale = F, verbose = T)
> sessionInfo()
R version 4.2.2 (2022-10-31)
Platform: x86_64-apple-darwin17.0 (64-bit)
Running under: macOS Ventura 13.3.1
Matrix products: default
LAPACK: /Library/Frameworks/R.framework/Versions/4.2/Resources/lib/libRlapack.dylib
locale:
[1] en_US.UTF-8/en_US.UTF-8/en_US.UTF-8/C/en_US.UTF-8/en_US.UTF-8
attached base packages:
[1] parallel stats4 stats graphics grDevices utils datasets methods base
other attached packages:
[1] doParallel_1.0.17 iterators_1.0.14 foreach_1.5.2 Pando_1.0.3 JASPAR2020_0.99.10
[6] TFBSTools_1.36.0 presto_1.0.0 data.table_1.14.8 harmony_0.1.1 Rcpp_1.0.10
[11] ensembldb_2.22.0 AnnotationFilter_1.22.0 GenomicFeatures_1.50.4 AnnotationDbi_1.60.2 Biobase_2.58.0
[16] simspec_0.0.0.9000 BSgenome.Hsapiens.UCSC.hg38_1.4.5 BSgenome_1.66.3 rtracklayer_1.58.0 Biostrings_2.66.0
[21] XVector_0.38.0 GenomicRanges_1.50.2 GenomeInfoDb_1.34.9 IRanges_2.32.0 S4Vectors_0.36.2
[26] AnnotationHub_3.6.0 BiocFileCache_2.6.1 dbplyr_2.3.2 BiocGenerics_0.44.0 Signac_1.9.0
[31] SeuratObject_4.1.3 Seurat_4.3.0 patchwork_1.1.2 lubridate_1.9.2 forcats_1.0.0
[36] stringr_1.5.0 dplyr_1.1.2 purrr_1.0.1 readr_2.1.4 tidyr_1.3.0
[41] tibble_3.2.1 ggplot2_3.4.2 tidyverse_2.0.0
loaded via a namespace (and not attached):
[1] Hmisc_5.1-0 ica_1.0-3 RcppRoll_0.3.0 ps_1.7.5 Rsamtools_2.14.0
[6] lmtest_0.9-40 rprojroot_2.0.3 crayon_1.5.2 MASS_7.3-60 nlme_3.1-162
[11] backports_1.4.1 qlcMatrix_0.9.7 rlang_1.1.1 ROCR_1.0-11 irlba_2.3.5.1
[16] callr_3.7.3 limma_3.54.2 filelock_1.0.2 BiocParallel_1.32.6 rjson_0.2.21
[21] CNEr_1.34.0 bit64_4.0.5 glue_1.6.2 poweRlaw_0.70.6 sctransform_0.3.5
[26] processx_3.8.1 vipor_0.4.5 chromVAR_1.20.2 spatstat.sparse_3.0-1 spatstat.geom_3.2-1
[31] tidyselect_1.2.0 SummarizedExperiment_1.28.0 usethis_2.1.6 fitdistrplus_1.1-11 XML_3.99-0.14
[36] zoo_1.8-12 GenomicAlignments_1.34.1 xtable_1.8-4 magrittr_2.0.3 evaluate_0.21
[41] cli_3.6.1 zlibbioc_1.44.0 rstudioapi_0.14 miniUI_0.1.1.1 sp_1.6-0
[46] rpart_4.1.19 fastmatch_1.1-3 maps_3.4.1 shiny_1.7.4 xfun_0.39
[51] pkgbuild_1.4.0 cluster_2.1.4 caTools_1.18.2 tidygraph_1.2.3 KEGGREST_1.38.0
[56] interactiveDisplayBase_1.36.0 ggrepel_0.9.3 biovizBase_1.46.0 listenv_0.9.0 TFMPvalue_0.0.9
[61] png_0.1-8 future_1.32.0 withr_2.5.0 bitops_1.0-7 slam_0.1-50
[66] ggforce_0.4.1 plyr_1.8.8 sparsesvd_0.2-2 pracma_2.4.2 pillar_1.9.0
[71] cachem_1.0.8 fs_1.6.2 hdf5r_1.3.8 vctrs_0.6.2 ellipsis_0.3.2
[76] generics_0.1.3 devtools_2.4.5 tools_4.2.2 foreign_0.8-84 beeswarm_0.4.0
[81] munsell_0.5.0 tweenr_2.0.2 DelayedArray_0.24.0 fastmap_1.1.1 compiler_4.2.2
[86] pkgload_1.3.2 abind_1.4-5 httpuv_1.6.11 sessioninfo_1.2.2 plotly_4.10.1
[91] GenomeInfoDbData_1.2.9 gridExtra_2.3 lattice_0.21-8 ggpointdensity_0.1.0 deldir_1.0-9
[96] utf8_1.2.3 later_1.3.1 jsonlite_1.8.4 scales_1.2.1 docopt_0.7.1
[101] pbapply_1.7-0 sparseMatrixStats_1.10.0 lazyeval_0.2.2 nabor_0.5.0 promises_1.2.0.1
[106] R.utils_2.12.2 goftest_1.2-3 spatstat.utils_3.0-3 reticulate_1.28 checkmate_2.2.0
[111] rmarkdown_2.21 cowplot_1.1.1 Rtsne_0.16 dichromat_2.0-0.1 uwot_0.1.14
[116] igraph_1.4.3 survival_3.5-5 yaml_2.3.7 htmltools_0.5.5 memoise_2.0.1
[121] VariantAnnotation_1.44.1 profvis_0.3.8 BiocIO_1.8.0 graphlayouts_1.0.0 viridisLite_0.4.2
[126] digest_0.6.31 mime_0.12 rappdirs_0.3.3 pals_1.7 RSQLite_2.3.1
[131] future.apply_1.11.0 mapproj_1.2.11 remotes_2.4.2 urlchecker_1.0.1 blob_1.2.4
[136] R.oo_1.25.0 splines_4.2.2 Formula_1.2-5 labeling_0.4.2 ProtGenerics_1.30.0
[141] RCurl_1.98-1.12 hms_1.1.3 colorspace_2.1-0 base64enc_0.1-3 BiocManager_1.30.20
[146] ggbeeswarm_0.7.2 ggrastr_1.0.1 nnet_7.3-19 RANN_2.6.1 ggseqlogo_0.1
[151] fansi_1.0.4 tzdb_0.4.0 parallelly_1.35.0 R6_2.5.1 grid_4.2.2
[156] ggridges_0.5.4 lifecycle_1.0.3 curl_5.0.0 leiden_0.4.3 motifmatchr_1.20.0
[161] Matrix_1.5-4.1 desc_1.4.2 RcppAnnoy_0.0.20 RColorBrewer_1.1-3 spatstat.explore_3.2-1
[166] htmlwidgets_1.6.2 polyclip_1.10-4 biomaRt_2.54.1 timechange_0.2.0 seqLogo_1.64.0
[171] globals_0.16.2 htmlTable_2.4.1 spatstat.random_3.1-5 progressr_0.13.0 codetools_0.2-19
[176] matrixStats_0.63.0 GO.db_3.16.0 gtools_3.9.4 prettyunits_1.1.1 R.methodsS3_1.8.2
[181] gtable_0.3.3 DBI_1.1.3 tensor_1.5 httr_1.4.6 KernSmooth_2.23-21
[186] stringi_1.7.12 progress_1.2.2 reshape2_1.4.4 farver_2.1.1 annotate_1.76.0
[191] viridis_0.6.3 DT_0.28 xml2_1.3.4 restfulr_0.0.15 scattermore_1.1
[196] BiocVersion_3.16.0 bit_4.0.5 MatrixGenerics_1.10.0 spatstat.data_3.0-1 ggraph_2.1.0
[201] pkgconfig_2.0.3 DirichletMultinomial_1.40.0 knitr_1.43
Not part of dev team but might be useful:
In my understanding infer_grn() will by default fit models for genes found as VariableFeatures. If VariableFeatures are not present in the object, it will throw the error that @TerezaClarence reported.
One can solve the problem by specifying genes using genes argument in the infer_grn().
