search

jackbibby1 / SCPA

R package for pathway analysis in scRNA-seq data
https://jackbibby1.github.io/SCPA/
GNU General Public License v3.0
57 stars 6 forks source link

Error reproducing step in "Systems level tissue comparison" tutorial #80

Closed JonasDBM closed 1 month ago

JonasDBM commented 1 month ago

Hi,

many thanks for developing and maintaining this very interesting tool! I tried to analyse my own data using the tutorial "Systems level tissue comparison" but ran into an error and got the same error when running the tutorial with the code and dataset provided.

Many thanks for your help!

Best,

Jonas

The specific error you're running into scpa_results <- Reduce(full_join, c(get_qvals(bl_bm, "bm"),

To Reproduce The code to reproduce the error e.g.

 library(SCPA)
  library(Seurat)
  library(msigdbr)
  library(magrittr)
  library(tidyverse)
  library(ComplexHeatmap)
pathways <- msigdbr(species = "Homo sapiens", category = "H") %>%
  format_pathways()
tissue_data <- szabo_t_cell
DimPlot(tissue_data, split.by = "tis_stim", group.by = "fine", ncol = 4)
cell_types <- unique(tissue_data$fine)
split_tissue <- SplitObject(tissue_data, split.by = "tissue")
# We first create empty lists to store results from the for loop
bl_bm <- list(); bl_ln <- list(); bl_lung <- list()
for (i in cell_types) {

  # We extract expression data using `seurat_extract` based on tissue, cell_type ("fine"), and stimulation ("none")
  blood <- seurat_extract(split_tissue$bl, 
                          meta1 = "fine", value_meta1 = i,
                          meta2 = "stimulation", value_meta2 = "none")

  bm <- seurat_extract(split_tissue$bm, 
                       meta1 = "fine", value_meta1 = i,
                       meta2 = "stimulation", value_meta2 = "none")

  ln <- seurat_extract(split_tissue$ln, 
                       meta1 = "fine", value_meta1 = i,
                       meta2 = "stimulation", value_meta2 = "none")

  lung <- seurat_extract(split_tissue$lung, 
                         meta1 = "fine", value_meta1 = i,
                         meta2 = "stimulation", value_meta2 = "none")

  # We then compare all tissues to the blood using `compare_pathways`
  print(paste("comparing", i))
  bl_bm[[i]] <- compare_pathways(list(blood, bm), pathways)
  bl_ln[[i]] <- compare_pathways(list(blood, ln), pathways)
  bl_lung[[i]] <- compare_pathways(list(blood, lung), pathways)

  # For faster analysis with parallel processing, use 'parallel = TRUE' and 'cores = x' arguments

}
get_qvals <- function(scpa_out, name) {

  df <- list()
  for (i in names(scpa_out)) {
    df[[i]] <- scpa_out[[i]] %>%
      select(Pathway, qval)
  }

  col_names <- names(df)
  for (i in 1:length(df)) {
    df[[i]] <- set_colnames(df[[i]], c("pathway", paste(name, col_names[[i]], sep = "_")))
  }

  return(df)

}
scpa_results <- Reduce(full_join, c(get_qvals(bl_bm, "bm"),
                                    get_qvals(bl_ln, "ln"),
                                    get_qvals(bl_lung, "lung")))

Additional context sessionInfo() R version 4.3.2 (2023-10-31) Platform: aarch64-apple-darwin20 (64-bit) Running under: macOS Ventura 13.2.1

Matrix products: default BLAS: /System/Library/Frameworks/Accelerate.framework/Versions/A/Frameworks/vecLib.framework/Versions/A/libBLAS.dylib LAPACK: /Library/Frameworks/R.framework/Versions/4.3-arm64/Resources/lib/libRlapack.dylib; LAPACK version 3.11.0

locale: [1] en_US.UTF-8/en_US.UTF-8/en_US.UTF-8/C/en_US.UTF-8/en_US.UTF-8

time zone: Europe/Zurich tzcode source: internal

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

other attached packages: [1] ComplexHeatmap_2.18.0 magrittr_2.0.3 dyno_0.1.2 dynwrap_1.2.4
[5] dynplot_1.1.2 dynmethods_1.0.5.9000 dynguidelines_1.0.1 dynfeature_1.0.1
[9] devtools_2.4.5 usethis_2.2.3 lubridate_1.9.3 forcats_1.0.0
[13] stringr_1.5.1 purrr_1.0.2 readr_2.1.5 tidyverse_2.0.0
[17] ggraph_2.2.1 igraph_2.0.3 org.Mm.eg.db_3.17.0 fastTopics_0.6-192
[21] Nebulosa_1.10.0 RColorBrewer_1.1-3 tradeSeq_1.14.0 condiments_1.8.0
[25] TSCAN_2.0.0 Lamian_0.99.2 slingshot_2.8.0 TrajectoryUtils_1.8.0
[29] princurve_2.1.6 viridis_0.6.5 viridisLite_0.4.2 annotables_0.2.0
[33] tibble_3.2.1 AnnotationDbi_1.62.2 scMiko_0.1.0 flexdashboard_0.6.2
[37] tidyr_1.3.1 gsdensity_0.1.3 msigdbr_7.5.1 SCPA_1.6.1
[41] SCpubr_2.0.2 singleCellTK_2.14.0 DelayedArray_0.28.0 SparseArray_1.2.3
[45] S4Arrays_1.2.0 abind_1.4-5 Matrix_1.6-5 patchwork_1.2.0
[49] scater_1.30.1 ggplot2_3.5.1 scuttle_1.12.0 miloR_1.8.1
[53] edgeR_4.0.11 limma_3.58.1 SingleCellExperiment_1.24.0 SummarizedExperiment_1.32.0 [57] Biobase_2.62.0 GenomicRanges_1.54.1 GenomeInfoDb_1.38.5 IRanges_2.36.0
[61] S4Vectors_0.40.2 BiocGenerics_0.46.0 MatrixGenerics_1.14.0 matrixStats_1.3.0
[65] dplyr_1.1.4 biomaRt_2.56.1 Seurat_5.1.0 lemur_1.3.2
[69] SeuratObject_5.0.2 sp_2.1-4

loaded via a namespace (and not attached): [1] pheatmap_1.0.12 nnls_1.5 DBI_1.2.3 bslib_0.7.0 httr_1.4.7
[6] ggh4x_0.2.8 BiocParallel_1.36.0 prettyunits_1.2.0 yulab.utils_0.1.4 ggplotify_0.1.2
[11] GenomicAlignments_1.36.0 sparseMatrixStats_1.14.0 brio_1.1.5 spatstat.geom_3.3-2 babelgene_22.9
[16] pillar_1.9.0 Rgraphviz_2.44.0 R6_2.5.1 mime_0.12 glmGamPoi_1.17.4
[21] reticulate_1.37.0 uwot_0.1.16 Rhdf5lib_1.24.1 ROCR_1.0-11 Hmisc_5.1-3
[26] downloader_0.4 parallelly_1.37.1 GlobalOptions_0.1.2 sparsesvd_0.2-2 caTools_1.18.2
[31] mgcv_1.9-1 polyclip_1.10-7 beachmat_2.18.0 htmltools_0.5.8.1 babelwhale_1.2.0
[36] fansi_1.0.6 caret_6.0-94 dnet_1.1.7 e1071_1.7-14 remotes_2.5.0
[41] ggrepel_0.9.5 fastICA_1.2-4 fgsea_1.26.0 spatstat.utils_3.0-5 HDO.db_0.99.1
[46] clusterProfiler_4.8.3 rpart_4.1.23 clue_0.3-65 sceasy_0.0.7 qlcMatrix_0.9.8
[51] scatterpie_0.2.3 fitdistrplus_1.2-1 goftest_1.2-3 tidyselect_1.2.1 RSQLite_2.3.7
[56] cowplot_1.1.3 GenomeInfoDbData_1.2.11 utf8_1.2.4 ScaledMatrix_1.10.0 scattermore_1.2
[61] sessioninfo_1.2.2 spatstat.data_3.1-2 gridExtra_2.3 fs_1.6.4 sctransform_0.4.1
[66] future.apply_1.11.2 graph_1.78.0 R.oo_1.26.0 topGO_2.52.0 ipred_0.9-15
[71] RcppHNSW_0.6.0 vipor_0.4.7 doRNG_1.8.6 scigenex_1.4.11 Rtsne_0.17
[76] DelayedMatrixStats_1.24.0 lazyeval_0.2.2 sass_0.4.9 scales_1.3.0 munsell_0.5.1
[81] lava_1.8.0 treeio_1.24.3 R.utils_2.12.3 profvis_0.3.8 pracma_2.4.4
[86] bitops_1.0-7 labeling_0.4.3 R.methodsS3_1.8.2 KEGGREST_1.40.1 promises_1.3.0
[91] shape_1.4.6.1 rhdf5filters_1.14.1 wesanderson_0.3.7 zoo_1.8-12 DropletUtils_1.22.0
[96] locfit_1.5-9.10 RSpectra_0.16-2 assertthat_0.2.1 paletteer_1.6.0 tools_4.3.2
[101] ape_5.8 processx_3.8.4 shiny_1.8.1.1 BiocFileCache_2.10.1 rlang_1.1.4
[106] generics_0.1.3 BiocSingular_1.18.0 ggridges_0.5.6 evaluate_0.24.0 fastcluster_1.2.6
[111] tictoc_1.2.1 reshape2_1.4.4 sciCSR_0.3.2 infotheo_1.2.0.1 expm_0.999-9
[116] colorspace_2.1-0 ellipsis_0.3.2 data.table_1.15.4 withr_3.0.0 RCurl_1.98-1.16
[121] xtable_1.8-4 plyr_1.8.9 invgamma_1.1 aplot_0.2.3 multimode_1.5
[126] ks_1.14.2 mclust_6.1.1 httpuv_1.6.15 rmarkdown_2.27 GSVAdata_1.36.0
[131] philentropy_0.8.0 MASS_7.3-60 dqrng_0.4.1 docopt_0.7.1 deldir_2.0-4
[136] GO.db_3.17.0 rhdf5_2.46.1 mixsqp_0.3-54 supraHex_1.38.0 tensor_1.5
[141] vctrs_0.6.5 lifecycle_1.0.4 proxy_0.4-27 codetools_0.2-20 fastDummies_1.7.3
[146] DT_0.33 recipes_1.0.10 nlme_3.1-165 dyndimred_1.0.4 combinat_0.0-8
[151] future_1.33.2 progress_1.2.3 dbplyr_2.5.0 pkgload_1.4.0 jquerylib_0.1.4
[156] Rcpp_1.0.13 SQUAREM_2021.1 rstudioapi_0.16.0 stringi_1.8.4 hms_1.1.3
[161] pbapply_1.7-2 cachem_1.1.0 pROC_1.18.5 BiocManager_1.30.23 hdf5r_1.3.11
[166] tidytree_0.4.6 listenv_0.9.1 XVector_0.42.0 urlchecker_1.0.1 plotly_4.10.4
[171] ggtree_3.8.2 enrichplot_1.20.3 GetoptLong_1.0.5 pkgbuild_1.4.4 ggfun_0.1.5
[176] HDF5Array_1.30.0 lmds_0.1.0 carrier_0.1.1 Formula_1.2-5 truncnorm_1.0-9
[181] clValid_0.7 htmlwidgets_1.6.4 kernlab_0.9-32 GA_3.2.4 markovchain_0.9.5
[186] matrixcalc_1.0-6 class_7.3-22 memoise_2.0.1 crayon_1.5.3 gridGraphics_0.5-1
[191] rappdirs_0.3.3 nbpMatching_1.5.5 xml2_1.3.6 filelock_1.0.3 GOSemSim_2.26.1
[196] ashr_2.2-63 png_0.1-8 progressr_0.14.0 tzdb_0.4.0 fastmap_1.2.0
[201] GSEABase_1.62.0 transport_0.15-2 tidygraph_1.3.1 pkgconfig_2.0.3 cli_3.6.3
[206] beeswarm_0.4.0 DOSE_3.26.2 ggforce_0.4.2 ps_1.7.7 prodlim_2024.06.25
[211] nnet_7.3-19 MAYA_1.0.0 lmtest_0.9-40 RcppArmadillo_0.12.8.4.0 RcppAnnoy_0.0.22
[216] amap_0.8-19 slam_0.1-51 timechange_0.3.0 foreign_0.8-87 timeDate_4032.109
[221] splines_4.3.2 askpass_1.2.0 blob_1.2.4 annotate_1.78.0 XML_3.99-0.17
[226] globals_0.16.3 ggbeeswarm_0.7.2 knitr_1.48 dynparam_1.0.2 distinct_1.12.2
[231] ggstar_1.0.4 ica_1.0-3 spam_2.10-0 compiler_4.3.2 rjson_0.2.21
[236] WriteXLS_6.7.0 RcppParallel_5.1.8 bit_4.0.5 hexbin_1.28.3 leidenbase_0.1.30
[241] diptest_0.77-1 BiocNeighbors_1.20.2 glue_1.7.0 umap_0.2.10.0 digest_0.6.36
[246] quadprog_1.5-8 irlba_2.3.5.1 leiden_0.4.3.1 graphlayouts_1.1.1 eds_1.2.0
[251] foreach_1.5.2 spatstat.random_3.3-1 SparseM_1.84-2 zlibbioc_1.48.0 GeneTrajectory_1.0.0
[256] dotCall64_1.1-1 tweenr_2.0.3 lattice_0.22-6 ModelMetrics_1.2.2.2 proxyC_0.4.1
[261] CelliD_1.11.1 statmod_1.5.0 openssl_2.2.0 rsvd_1.0.5 gson_0.1.0
[266] yaml_2.3.9 mvtnorm_1.2-5 qvalue_2.32.0 later_1.3.2 backports_1.5.0
[271] shadowtext_0.1.4 Rsamtools_2.16.0 parallel_4.3.2 rematch2_2.1.2 miniUI_0.1.1.1
[276] gtable_0.3.5 xfun_0.46 Biostrings_2.68.1 curl_5.2.1 multicross_2.1.0
[281] rootSolve_1.8.2.4 doParallel_1.0.17 org.Hs.eg.db_3.17.0 KernSmooth_2.23-24 survival_3.7-0
[286] desc_1.4.3 jsonlite_1.8.8 harmony_1.2.0 base64enc_0.1-3 iterators_1.0.14
[291] spatstat.univar_3.0-0 crossmatch_1.4-0 RhpcBLASctl_0.23-42 testthat_3.2.1.1 checkmate_2.3.1
[296] fastmatch_1.1-4 gower_1.0.1 hardhat_1.3.1 iheatmapr_0.7.1 gtools_3.9.5
[301] htmlTable_2.4.3 spatstat.sparse_3.1-0 rngtools_1.5.2 RANN_2.6.1 circlize_0.4.16
[306] spatstat.explore_3.3-1 Ecume_0.9.2 anticlust_0.8.5 bit64_4.0.5 cluster_2.1.6
[311] farver_2.1.2 dynutils_1.0.11 ranger_0.16.0 gplots_3.1.3.1

jackbibby1 commented 1 month ago

Hey,

Thanks for the info. It looks like this may be an issue with masking of the select function (another library you have loaded also has a select function that's overriding the dplyr one). If you replace select with dplyr::select, hopefully that should be solved.

Thanks Jack

JonasDBM commented 1 month ago

Thanks for the fast help!

dplyr:: select and magrittr:set_colnames fixed it!

Best,

Jonas