library(ggplot2)
library(vctrs)
library(dplyr)
library(patchwork, lib.loc = "/apps/rocs/2020.08/cascadelake/software/R/4.1.2-foss-2020a/lib64/R/library")
library(Seurat)
library(dorothea)
library(tibble)
library(readr)
library(stringr)
library(cowplot)
library(RColorBrewer)
library(mistyR)
data_directory <- params$data_directory
analysis_name <- params$analysis_name
input_folder <- "Cell2Location/inputs/scRNAseq-ref/raw"
files_to_read <- paste0(data_directory, analysis_name,input_folder)
Reading Seurat objects
expression_matrix <- Read10X(data.dir = files_to_read)
seurat_object <-
CreateSeuratObject(counts = expression_matrix, min.cells = 10, min.features = 200)
metedata_df <- read_tsv(paste0(files_to_read, "/metadata.tsv")) %>%
tibble::column_to_rownames(var="CELL")
seurat_object <- AddMetaData(seurat_object, metedata_df)
seurat_object <- SCTransform(seurat_object, assay = "RNA", verbose = FALSE)
According to Lee’s publication, we select patients with a CMS2 phenotype.
CMS2_patients <- c("SMC22-T", "SMC18-T", "SMC21-T", "SMC09-T",
"SMC23-T", "SMC25-T", "SMC11-T", "SMC07-T")
seurat_object_CMS2 <- subset(x=seurat_object, idents = CMS2_patients,
subset = Cell_subtype != "Unknown")
Idents(seurat_object_CMS2) <- seurat_object_CMS2$Cell_subtype
seurat_object_CMS2 <- subset(seurat_object_CMS2)
We run dorothea with the same parameters as in the ST data.
dorothea_regulon <-
get(data("dorothea_hs", package = "dorothea"))
confidence_levels<- c("A","B","C")
regulon_filtered <- dorothea_regulon %>%
dplyr::filter(confidence %in% confidence_levels)
seurat_object_CMS2 <-
run_viper(seurat_object_CMS2, regulon_filtered, assay_key = "SCT",
options = list(method = "scale", minsize = 4, eset.filter = FALSE,
cores = 1, verbose = FALSE))
We load Misty results to get the TF of interest
results_misty_folders <- paste0(data_directory, analysis_name,
"IntermediaryFiles/Misty_Results/results_TF_Ligands_DorotheaClusters/")
misty_results <-
collect_results(as.list(list.dirs(results_misty_folders, recursive = FALSE)))
TFs_toconsider <-
misty_results$importances.aggregated %>%
dplyr::filter(view == "juxta.ligands_2") %>%
dplyr::filter(Importance > 1) %>%
dplyr::pull(Target) %>% unique()
We plot the average TF activity per cell type
DotPlot(seurat_object_CMS2, assay = "dorothea",
features = TFs_toconsider, dot.scale = 5,
col.min = 0) +
RotatedAxis() + labs(y= "Cell Types", x = "TFs") +
scale_color_gradient(low = "white", high = "darkred") +
guides(color = guide_colorbar(title = "Average TF activity"),
size = guide_legend("Percent Active"))
## R version 4.1.2 (2021-11-01)
## Platform: x86_64-pc-linux-gnu (64-bit)
## Running under: Fedora 33 (Container Image)
##
## Matrix products: default
## BLAS/LAPACK: /apps/rocs/2020.08/cascadelake/software/OpenBLAS/0.3.9-GCC-9.3.0/lib/libopenblas_skylakexp-r0.3.9.so
##
## locale:
## [1] LC_CTYPE=en_US.UTF-8 LC_NUMERIC=C
## [3] 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
## [7] LC_PAPER=en_US.UTF-8 LC_NAME=C
## [9] LC_ADDRESS=C LC_TELEPHONE=C
## [11] LC_MEASUREMENT=en_US.UTF-8 LC_IDENTIFICATION=C
##
## attached base packages:
## [1] stats graphics grDevices utils datasets methods base
##
## other attached packages:
## [1] purrr_0.3.4 mistyR_1.2.1 RColorBrewer_1.1-3
## [4] cowplot_1.1.1 stringr_1.4.0 readr_2.1.2
## [7] tibble_3.1.7 dorothea_1.6.0 sp_1.5-0
## [10] SeuratObject_4.1.0 Seurat_4.1.0 patchwork_1.1.1
## [13] dplyr_1.0.9 vctrs_0.4.1 ggplot2_3.3.6
## [16] BiocManager_1.30.18
##
## loaded via a namespace (and not attached):
## [1] plyr_1.8.7 igraph_1.3.2 lazyeval_0.2.2
## [4] splines_4.1.2 listenv_0.8.0 scattermore_0.8
## [7] digest_0.6.29 htmltools_0.5.2 fansi_1.0.3
## [10] magrittr_2.0.3 tensor_1.5 cluster_2.1.3
## [13] mixtools_1.2.0 ROCR_1.0-11 tzdb_0.3.0
## [16] globals_0.15.0 matrixStats_0.62.0 R.utils_2.11.0
## [19] vroom_1.5.7 spatstat.sparse_2.1-1 colorspace_2.1-0
## [22] ggrepel_0.9.1 xfun_0.31 crayon_1.5.1
## [25] jsonlite_1.8.0 progressr_0.10.1 spatstat.data_2.2-0
## [28] survival_3.3-1 zoo_1.8-10 glue_1.6.2
## [31] polyclip_1.10-0 gtable_0.3.0 leiden_0.4.2
## [34] kernlab_0.9-31 future.apply_1.9.0 BiocGenerics_0.40.0
## [37] abind_1.4-7 scales_1.2.0 DBI_1.1.3
## [40] spatstat.random_2.2-0 miniUI_0.1.1.1 Rcpp_1.0.8.3
## [43] viridisLite_0.4.0 xtable_1.8-6 reticulate_1.25
## [46] spatstat.core_2.4-4 bit_4.0.4 proxy_0.4-27
## [49] htmlwidgets_1.5.4 httr_1.4.3 ellipsis_0.3.2
## [52] ica_1.0-2 farver_2.1.0 R.methodsS3_1.8.2
## [55] pkgconfig_2.0.3 sass_0.4.1 uwot_0.1.11
## [58] deldir_1.0-6 utf8_1.2.2 labeling_0.4.2
## [61] tidyselect_1.1.2 rlang_1.0.2 reshape2_1.4.4
## [64] later_1.3.0 munsell_0.5.0 tools_4.1.2
## [67] cli_3.3.0 generics_0.1.2 ggridges_0.5.3
## [70] evaluate_0.15 fastmap_1.1.0 yaml_2.3.5
## [73] goftest_1.2-3 knitr_1.39 bit64_4.0.5
## [76] fitdistrplus_1.1-8 bcellViper_1.30.0 RANN_2.6.1
## [79] pbapply_1.5-0 future_1.26.1 nlme_3.1-158
## [82] mime_0.12 R.oo_1.25.0 compiler_4.1.2
## [85] rstudioapi_0.13 plotly_4.10.0 png_0.1-7
## [88] e1071_1.7-11 spatstat.utils_2.3-1 viper_1.28.0
## [91] bslib_0.3.1 stringi_1.7.6 highr_0.9
## [94] rgeos_0.5-10 lattice_0.20-45 Matrix_1.4-2
## [97] furrr_0.3.0 pillar_1.7.0 lifecycle_1.0.1
## [100] spatstat.geom_2.4-0 lmtest_0.9-40 jquerylib_0.1.4
## [103] RcppAnnoy_0.0.19 data.table_1.14.2 irlba_2.3.5
## [106] httpuv_1.6.5 R6_2.5.1 promises_1.2.0.1
## [109] KernSmooth_2.23-20 gridExtra_2.3 parallelly_1.32.0
## [112] codetools_0.2-18 MASS_7.3-57 assertthat_0.2.1
## [115] withr_2.5.0 sctransform_0.3.3 mgcv_1.8-40
## [118] parallel_4.1.2 hms_1.1.1 grid_4.1.2
## [121] rpart_4.1.16 tidyr_1.2.0 class_7.3-20
## [124] rmarkdown_2.14 segmented_1.6-0 Rtsne_0.16
## [127] Biobase_2.54.0 shiny_1.7.1