NB-FOXR2 project [source]
Configuration of project directory & analysis outputs:
Show full config
source(here("rr_helpers.R"))
# Set up outputs
message("Document index: ", doc_id)
## Document index: 06.1
# Specify where to save outputs
out <- here("output", doc_id); dir.create(out, recursive = TRUE)
figout <- here("figures", doc_id, "/")
cache <- paste0("/project/kleinman/bhavyaa.chandarana/cache/NB-FOXR2/public/", doc_id, "/")
The root directory of this project is:
## /project/kleinman/bhavyaa.chandarana/from_hydra/2023-05-NB-FOXR2/public
Outputs and figures will be saved at these paths, relative to project root:
## public/output/06.1
## public/figures/06.1//
Setting a random seed:
set.seed(100)
This document explores projections/mapping of tumors to normal using both bulk and single-cell RNAseq data.
With bulk RNAseq data, the strategy is to evaluate enrichment of cell type-specific gene signatures (derived in 02
) in each bulk transcriptome. We can then evaluate the highest scoring signature(s) in each sample, or compare enrichment scores of gene signatures between tumor types.
library(here)
library(magrittr)
library(tidyr)
library(dplyr)
library(readr)
library(readxl)
library(stringr)
library(glue)
library(purrr)
library(ggplot2)
library(ggrepel)
library(tibble)
library(broom)
library(cowplot)
library(Seurat)
library(ComplexHeatmap)
library(data.table)
library(ggExtra)
library(Rtsne)
source(here("include/style.R"))
source(here("code/functions/RNAseq.R"))
source(here("code/functions/scRNAseq.R"))
source(here("code/functions/ssGSEA.R"))
ggplot2::theme_set(theme_min())
square_theme <- theme(aspect.ratio = 1)
large_text <- theme(text = element_text(size = 15))
conflicted::conflicts_prefer(dplyr::filter)
# Function to load an RData file and return it
loadRData <- function(fileName){
load(fileName)
get(ls()[ls() != "fileName"])
}
meta <- read_tsv(here("output/00/metadata_patients_NGS.tsv"))
## Rows: 176 Columns: 20
## ── Column specification ────────────────────────────────────────────────────────
## Delimiter: "\t"
## chr (20): Sample, ID_Patient, ID_Sample, ID, FOXR2_positive, Group, Source, ...
##
## ℹ Use `spec()` to retrieve the full column specification for this data.
## ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.
meta_bulk <- meta %>% filter(RNAseq == "Y")
Load cell type specific gene signatures, compiled in document 02
, to use as input for ssGSEA:
load(here("output/02/signatures.Rda"))
Load the raw counts from the bulk RNAseq pipeline, for the test set of pediatric brain tumors:
# load counts aggregated with in-house pipeline
test_raw_counts <- file.path(here("data/RNAseq/pipeline_l3/2023-05-test_pbt/counts/Ensembl.ensGene.exon.raw.tsv.gz")) %>%
read.table(header = T, row.names = 1, check.names = F, sep = "\t")
# rename genes with ENSEMBL ID only
head(rownames(test_raw_counts) %<>% strsplit(":") %>% sapply(getElement, 1))
## [1] "ENSG00000118473" "ENSG00000162426" "ENSG00000157191" "ENSG00000169504"
## [5] "ENSG00000142920" "ENSG00000186094"
test_raw_counts <- as.matrix(test_raw_counts)
test_raw_counts[1:5, 1:5]
## 00-009P 00-011R 00-020P 00-027P 00-029P
## ENSG00000118473 118 1400 1873 1978 1274
## ENSG00000162426 16 345 482 839 558
## ENSG00000157191 947 2035 3668 6433 6711
## ENSG00000169504 1947 14701 31704 19765 29772
## ENSG00000142920 142 394 930 1465 1241
Load raw counts for the background set of brain tumors/normal brain:
# load counts aggregated with in-house pipeline
background_raw_counts <- file.path(here("data/RNAseq/pipeline_l3/2023-04-background_samples/counts/Ensembl.ensGene.exon.raw.tsv.gz")) %>%
read.table(header = T, row.names = 1, check.names = F, sep = "\t")
# rename genes with ENSEMBL ID only
head(rownames(background_raw_counts) %<>% strsplit(":") %>% sapply(getElement, 1))
## [1] "ENSG00000118473" "ENSG00000162426" "ENSG00000157191" "ENSG00000169504"
## [5] "ENSG00000142920" "ENSG00000186094"
background_raw_counts <- as.matrix(background_raw_counts)
background_raw_counts[1:5, 1:5]
## EBT_1 EBT_2 EBT_3 EBT_4 EBT_5
## ENSG00000118473 384 557 54 160 515
## ENSG00000162426 154 128 199 176 118
## ENSG00000157191 2205 1570 1133 992 1065
## ENSG00000169504 7780 3180 9300 5670 4520
## ENSG00000142920 303 246 512 230 182
For extra-cranial neuroblastomas from Gartlguber et al, Nature Cancer, 2021 (n=579) processed data was obtained from the Shiny App associated with the study.
# load raw counts provided by the authors
ecnb_counts <- readRDS(here("data/RNAseq/external_data/Gartlgruber_NatCancer_2021/tumor_RNAseq_Counts_Matrix.rds"))
head(rownames(ecnb_counts))
## [1] "ENSG00000223972.4|DDX11L1" "ENSG00000227232.4|WASH7P"
## [3] "ENSG00000243485.2|MIR1302-11" "ENSG00000237613.2|FAM138A"
## [5] "ENSG00000268020.2|OR4G4P" "ENSG00000240361.1|OR4G11P"
# extract ENSEMBL id only from gene names
head(rownames(ecnb_counts) %<>% strsplit("\\.") %>% sapply(getElement, 1))
## [1] "ENSG00000223972" "ENSG00000227232" "ENSG00000243485" "ENSG00000237613"
## [5] "ENSG00000268020" "ENSG00000240361"
ecnb_counts <- as.matrix(ecnb_counts)
ecnb_counts[1:5, 1:5]
## NSP001-PT01 NSP002-PT01 NSP003-RM01 NSP005-RT01 NSP006-RM01
## ENSG00000223972 4 0 0 3 1
## ENSG00000227232 431 1348 1242 1117 1391
## ENSG00000243485 0 0 0 0 0
## ENSG00000237613 0 0 0 0 0
## ENSG00000268020 0 0 0 0 0
Load the sample metadata for Gartlgruber EC-NB dataset:
meta_ecnb <- readRDS(here("data/RNAseq/external_data/Gartlgruber_NatCancer_2021/annotation_tumor_RNAseq.rds"))
# load annotation of ChIPseq data to look for overlapping samples
meta_ecnb_chip <- readRDS(here("data/ChIPseq/external_data/Gartlgruber_NatCancer_2021/annotation_tumor_ChIPseq.rds"))
meta_ecnb <- meta_ecnb %>%
left_join(meta_ecnb_chip, by = c("ProjectID", "MYCN", "Stage", "Age", "Risk", "Relapse")) %>%
dplyr::rename(ID = ProjectID) %>%
mutate(Group = case_when(
MYCN == "Amp" ~ "EC-NB (MYCN-amp)",
TRUE ~ "EC-NB"))
Add metadata column for FOXR2 expression status in EC-NB:
# Load files and remove unused ones to save space
load(here("output/03/Gartlgruber_et_al_counts.Rda"))
rm(ecnb_counts_norm)
rm(ecnb_counts_vst)
# DESeq norm expression threshold used for FOXR2 +/-
threshold <- 2
foxr2_pos_samples <- ecnb_counts_tidy %>%
filter(gene_symbol == "FOXR2") %>%
filter(gene_expression > threshold) %>%
.$sample
meta_ecnb <- meta_ecnb %>% mutate(FOXR2_positive = case_when(ID %in% foxr2_pos_samples ~ "Y",
T ~ "N")) %>%
mutate(Group = "EC-NB") %>%
mutate(Source = "Gartlgruber et al. 2021")
rm(ecnb_counts_tidy)
Combine with brain tumor dataset:
meta_all <- bind_rows(meta_bulk, meta_ecnb)
Run ssGSEA for each dataset:
# test set
bulk_ssgsea <- compute_scores_bin(expr_mat = test_raw_counts,
gene_sets = signatures_ens_all,
save_le = FALSE,
alpha = 0.75,
normalize = FALSE,
n_cores = 1)
saveRDS(bulk_ssgsea, file = glue("{out}/ssgsea_test.Rds"))
# background set
bulk_ssgsea_bg <- compute_scores_bin(expr_mat = background_raw_counts,
gene_sets = signatures_ens_all,
save_le = FALSE,
alpha = 0.75,
normalize = FALSE,
n_cores = 1)
saveRDS(bulk_ssgsea_bg, file = glue("{out}/ssgsea_background.Rds"))
# extra-cranial NB
bulk_ssgsea_ecnb <- compute_scores_bin(expr_mat = ecnb_counts,
gene_sets = signatures_ens_all,
save_le = FALSE,
alpha = 0.75,
normalize = FALSE,
n_cores = 1)
saveRDS(bulk_ssgsea_ecnb, file = glue("{out}/ssgsea_ecnb.Rds"))
First, we will evaluate the ssGSEA scores for reference signatures in our background set of tumors and normal brain. To do this, we count the number of samples where each signature appears in the top 10. We'll consider the top decile as indiscriminate.
# control set
bulk_ssgsea_bg <- readRDS(glue("{out}/ssgsea_background.Rds"))
dim(bulk_ssgsea_bg)
## [1] 412 211
# filter to retained signatures
bulk_ssgsea_bg_keep <- bulk_ssgsea_bg %>%
filter(Signature %in% signatures_keep)
dim(bulk_ssgsea_bg_keep)
## [1] 390 211
# tidy dataframe
bulk_ssgsea_bg_top <- bulk_ssgsea_bg_keep %>%
gather(Sample, Score, 2:ncol(.)) %>%
group_by(Sample) %>%
top_n(10, Score)
# get n signatures - this many signatures appear in the top 10 of any sample
n_sigs_top <- length(unique(bulk_ssgsea_bg_top$Signature))
n_sigs_top
## [1] 135
# count the number of times each signature appears in the top 10 signatures for any sample
counts_top <- bulk_ssgsea_bg_top %>%
group_by(Signature) %>%
dplyr::count()
# show the values
counts_top %>% arrange(desc(n)) %>% DT::datatable()
# get threshold at top decile
(q <- quantile(counts_top$n, 0.9)[[1]])
## [1] 41
Then, we can identify indiscriminant signatures by calculating a threshold \(q\) where the probability that \(X\) is less than or equal to \(q\) is 90%.
counts_top %>%
ggplot(aes(x = n)) +
geom_histogram() +
geom_vline(xintercept = q, color = "red") +
xlab("Number of samples where the signatures is in the top 10") +
ylab("Number of signatures")
## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.
In other words, we identify indiscriminant signatures as those which are in the top 10 of more samples than 90% of other signatures. As a sanity check, we have 135 signatures appearing in the top 10 of any sample, and then are taking the top 13.5 signatures as indiscriminant.
# indiscriminant signatures
(sigs_indiscriminant <- counts_top %>% filter(n >= q) %>% pull(Signature))
## [1] "CGE-ANK3+,ENC1+" "CGE-GRIA1+,SST+"
## [3] "F-e15_CEXN2" "HF OPC"
## [5] "HP AST-FB" "HP Neu-NRGN-I"
## [7] "HP Neu-NRGN-II" "HP OPC"
## [9] "HP Oligodendrocytes" "Human fetal GE OPC"
## [11] "Human fetal GE P3" "Human fetal cerebellum 09-Ast"
## [13] "Str/SVZ 33-Neurogenic progenitor" "VSVZ-1 Neuron"
## [15] "VSVZ-2 Neuron"
length(sigs_indiscriminant)
## [1] 15
In addition, we also filter out a MGE signature which we determined to be not specific to its associated cell type, based on expression of the genes in this signature across multiple cell types. (See following section.)
sigs_indiscriminant <- c("MGE-NR2F1+,MEIS2+", sigs_indiscriminant)
length(sigs_indiscriminant)
## [1] 16
save(sigs_indiscriminant, file = glue("{out}/sigs_indiscriminant.Rda"))
In previous iterations of this analysis, MGE signature MGE-NR2F1+,MEIS2+
appeared in top ssGSEA matches in the majority of DIPG-H3K27M. This was unexpected because DIPG-H3K27M appear in the pons, not in the forebrain, where the MGE structure appears during development.
Let us check the detection rate of each gene of this signature in cells from a normal developmental atlas, namely MGE neurons vs. cell types that give rise to DIPG (OPC and oligodendrocytes) in both the pons and forebrain. We will also compare this with non-neural cell types (vascular, immune, endothelial) as a negative control.
Define a function to plot detection rate per gene in a given signature:
# Load full cell type ontology table for mouse developing forebrain and pons
ontology_jessa_full <- fread("/project/kleinman/selin.jessa/from_narval/HGG-oncohistones/public/output/05/TABLE_mouse_cluster_info.tsv") %>%
select(Label, Level3_type) %>%
rename(Label_broad = Level3_type)
# Function producing box plot of detection rate per gene in a signature
# Within cell populations: MGE neurons, Forebrain OPC & oligodendrocytes, Pons OPC & oligodendrocytes, and non-neural cells
plot_mge_neur_opc_box_plot <- function(signature,
signature_id,
conv_hg_to_mm){
print(signature_id)
ontology <- ontology_jessa_full
if(conv_hg_to_mm){
# Convert genes from human to mouse using cached biomaRt reference
# Produces file "genes_lds"
load(here("data/singlecell/references_genome/biomaRt_mm_to_hg_lds.Rda"))
# Store converted signature
mm_signature <- genes_lds %>% filter(HGNC.symbol %in% signature)
# Print any human genes that could not be converted
if(length(mm_signature) != length(signature)){
not_conv <- setdiff(signature, mm_signature$HGNC.symbol)
print(glue("The following genes of the provided signature could not be converted from human to mouse using biomaRt:"))
print(list(not_conv))
}
# Replace user-provided human signature with signature in mouse
signature <- mm_signature$MGI.symbol
rm(genes_lds)
}
# Load developing mouse atlas objects
atlas_dir <- "/project/kleinman/selin.jessa/from_hydra/atlas/data"
seurat_ct <- loadRData(glue("{atlas_dir}/joint_cortex_extended/joint_cortex_extended.seurat_v3.Rda"))
suppressMessages({
seurat_ct <- UpdateSeuratObject(seurat_ct)
})
seurat_po <- loadRData(glue("{atlas_dir}/joint_pons_extended/joint_pons_extended.seurat_v3.Rda"))
suppressMessages({
seurat_po <- UpdateSeuratObject(seurat_po)
})
# Store forebrain and pons objects in list and remove the originals
seurat_list <- c(seurat_ct, seurat_po)
names(seurat_list) <- c("cortex", "pons")
rm(seurat_ct)
rm(seurat_po)
# Add broad labels and categories for plot to Seurat objects
for(seurat_id in names(seurat_list)){
seurat <- seurat_list[[seurat_id]]
# Get metadata from Seurat object
labels <- seurat@meta.data %>%
select(ID_20210710) %>%
rownames_to_column(var = "cellname")
rownames(labels) <- NULL
# Add broad labels
labels <- labels %>%
left_join(ontology %>% select(Label, Label_broad),
by = c('ID_20210710' = 'Label')) %>%
column_to_rownames(var = "cellname")
seurat <- AddMetaData(seurat, labels)
# Group cells of selected cell type(s) for different rows of plot
if(seurat_id == "cortex"){
labels_ct <- labels %>%
mutate(plot_row = case_when(Label_broad %in% c("MGE inhibitory neurons") ~ "MGE inhibitory",
Label_broad %in% c("OPC", "Oligodendrocytes") ~ "OPC/Oligo (cortex)",
Label_broad %in% c("Meninges", "Endothelial", "Immune", "Pericytes") ~ "Endothelial/Immune (cortex)",
TRUE ~ "Other"))
seurat <- AddMetaData(seurat, labels_ct)
} else if(seurat_id == "pons"){
labels_po <- labels %>%
mutate(plot_row = case_when(Label_broad %in% c("OPC", "Oligodendrocytes") ~ "OPC/Oligo (pons)",
TRUE ~ "Other"))
seurat <- AddMetaData(seurat, labels_po)
}
seurat_list[[seurat_id]] <- seurat
}
# Retrieve detection rate for each gene in the signature
seurat_list[["cortex"]] <- seurat_list[["cortex"]][signature, ]
seurat_list[["pons"]] <- seurat_list[["pons"]][signature, ]
mat_list <- c()
mat_list[["cortex"]] <- as.matrix(seurat_list[["cortex"]]@assays$RNA@data)
mat_list[["pons"]] <- as.matrix(seurat_list[["pons"]]@assays$RNA@data)
pct_list <- c()
for(name in names(mat_list)){
x <- mat_list[[name]]
# Binarize
x[x > 0] <- 1
# Transpose preserving cell names
cellname <- colnames(x)
x <- as.data.table(t(x))
rownames(x) <- cellname
# Add in_group info for cells
x[, plot_row := as.character(seurat_list[[name]]@meta.data[["plot_row"]])]
# Get prop of cells expressing a gene, within group and out of group
pct <- x[, lapply(.SD, prop), by = plot_row]
# Convert to tidy format for plotting
pct <- pct %>%
pivot_longer(cols = colnames(pct)[-1], names_to = "Gene")
pct_list[[name]] <- pct
}
# Clean up
rm(seurat_list)
rm(mat_list)
# Join cortex and pons detection rates, removing "Other" cell type category from each
pct_list[["cortex"]] <- pct_list[["cortex"]] %>%
as.data.frame() %>%
filter(plot_row != "Other") %>%
rename(Cell_class = plot_row)
pct_list[["pons"]] <- pct_list[["pons"]] %>%
as.data.frame() %>%
filter(plot_row != "Other") %>%
rename(Cell_class = plot_row)
pct <- rbind(pct_list[["cortex"]], pct_list[["pons"]])
suppressMessages({
# Produce horizontal box plot
pct %>%
mutate(Cell_class = factor(Cell_class,
levels = rev(c("MGE inhibitory",
"OPC/Oligo (cortex)",
"OPC/Oligo (pons)",
"Endothelial/Immune (cortex)")))) %>%
ggplot(aes(x = value, y = Cell_class)) +
ggtitle(glue("{signature_id}")) +
geom_boxplot() +
xlab("% of cells expressing gene") +
ylab("") +
scale_x_continuous(breaks = c(0, 0.25, 0.5, 0.75, 1),
limits = c(-0.05, 1.05)) +
theme(axis.title.x = element_text(hjust = 0.5))
})
}
Plot detection rate in the signature of interest MGE-NR2F1+,MEIS2+
, along with two other MGE signatures used in the study:
sigs <- c("MGE-NR2F1+,MEIS2+",
"F-e12_MGINH",
"Human fetal GE MGE1")
plot_list <- map(sigs,
~ plot_mge_neur_opc_box_plot(signature = signatures_sym_all[[.x]],
signature_id = .x,
conv_hg_to_mm = T) + large_text)
## [1] "MGE-NR2F1+,MEIS2+"
## The following genes of the provided signature could not be converted from human to mouse using biomaRt:
## [[1]]
## [1] "NNAT" "BEX3" "HMGN2" "TMSB15A" "PPIA" "VDAC1" "TCEAL2"
## [8] "LDHA" "RACK1"
##
## [1] "F-e12_MGINH"
## The following genes of the provided signature could not be converted from human to mouse using biomaRt:
## [[1]]
## [1] "NNAT" "H3F3B" "HIST3H2BB" "ATPIF1" "AES" "USMG5"
## [7] "AHSA2" "UBL5" "DNAJA1" "MAP1LC3B2"
##
## [1] "Human fetal GE MGE1"
## The following genes of the provided signature could not be converted from human to mouse using biomaRt:
## [[1]]
## [1] "DLX6-AS1" "CADPS" "CEMIP2" "RIPOR2" "TMEM161B-AS1"
## [6] "C1orf61" "MIR181A1HG" "AES" "HIST1H2BG" "MALAT1"
## [11] "MEG3" "AKR1C1" "KCNQ1OT1"
plot_grid(plotlist = plot_list, nrow = 1)
The MGE-NR2F1+,MEIS2+
gene signature actually has a lower detection rate in MGE inhibitory neurons than the OPC of the forebrain and pons. Therefore, we remove this signature from our analysis for non-specific expression to the MGE.
Here, we'll prepare a supplementary table containing all the reference gene signatures and their associated cell types/datasets from which they originate.
# collapse signatures into strings
sigs_sym_df <- imap_dfr(signatures_sym_all,
~ data.frame(Cluster = .y, Signature_sym = glue_collapse(.x, sep = ",")))
sigs_ens_df <- imap_dfr(signatures_ens_all,
~ data.frame(Cluster = .y, Signature_ENS = glue_collapse(.x, sep = ",")))
# join with cell type annotation
TABLE_ref_sigs <- cell_type_anno_all %>%
left_join(sigs_sym_df) %>%
left_join(sigs_ens_df) %>%
# filter to kept signatures
filter(Keep) %>%
select(-Keep) %>%
mutate(Cell_type = ifelse(is.na(Cell_type), Cluster, Cell_type)) %>%
mutate(Included = ifelse(Cluster %in% sigs_indiscriminant, "N", "Y"))
## Joining with `by = join_by(Cluster)`
## Joining with `by = join_by(Cluster)`
rr_write_tsv(TABLE_ref_sigs,
glue("{out}/TABLE_reference_signatures.tsv"),
"Cell type specific reference signatures and associated metadata and dataset info")
## ...writing description of TABLE_reference_signatures.tsv to public/output/06.1/TABLE_reference_signatures.desc
Combine the ssGSEA scores from the two datasets:
bulk_ssgsea_test <- readRDS(glue("{out}/ssgsea_test.Rds"))
bulk_ssgsea_ecnb <- readRDS(glue("{out}/ssgsea_ecnb.Rds"))
scores_all <- full_join(bulk_ssgsea_test, bulk_ssgsea_ecnb, by = "Signature") %>%
as.data.frame() %>%
# filter out short signatures
filter(Signature %in% signatures_keep) %>%
# filter out the indiscriminant signature
filter(!(Signature %in% sigs_indiscriminant))
scores_all$Signature %>% unique %>% length()
## [1] 374
# correct a sample naming issue, where samples with numeric names XX get converted to XX.0
# and save a map between the coerced and cleaned names
sample_naming_map <- data.frame(Sample = colnames(scores_all))
sample_idx <- which(grepl("\\.0$", colnames(scores_all)))
colnames(scores_all)[sample_idx] <- gsub("\\.0$", "", colnames(scores_all)[sample_idx])
sample_naming_map$Sample_clean <- colnames(scores_all)
# get the set of test samples
test_samples <- base::setdiff(colnames(bulk_ssgsea_test), "Signature")
test_samples <- plyr::mapvalues(test_samples, from = sample_naming_map$Sample, to = sample_naming_map$Sample_clean, warn_missing = FALSE)
scores_tidy <- scores_all %>%
# convert to tidy/long format
gather("Sample", "Score", 2:ncol(.)) %>%
# join with sample metadata
left_join(meta_all, by = c("Sample" = "ID")) %>%
# join with signature metadata
left_join(cell_type_anno_all %>% dplyr::rename("Reference_sample" = Sample, "Reference_age" = Age),
by = c("Signature" = "Cluster")) %>%
# set ordering
mutate(Group = factor(Group, levels = names(palette_groups))) %>%
arrange(Group) %>%
mutate(Sample = factor(Sample, levels = unique(.$Sample))) %>%
mutate(Class = factor(Class, levels = names(palette_class)))
save(scores_all, scores_tidy, test_samples, file = glue("{out}/ssgsea_scores_tidy.Rda"))
# count samples
scores_tidy %>% distinct(Sample, Group) %>% group_by(Group) %>% count()
In-house processed data for PBT alone:
perp <- 10
tsne_plot <- scores_all[, c("Signature", test_samples)] %>%
tibble::column_to_rownames(var = "Signature") %>%
as.matrix() %>%
ssgsea_dim_red(method = "tsne", return_df = TRUE, perplexity = perp) %>%
left_join(meta_all, by = c("Sample" = "ID")) %>%
plot_dim_red(method = "tsne", palette = palette_groups,
point_size = 3, alpha = 0.8, title = glue("Perplexity: {perp}")) + square_theme + large_text
## @ returning data...
tsne_plot
PBT with Gartlgruber dataset High risk & Stage 4 EC-NB. (EC-NB samples with NA in either Risk or Stage metadata columns are removed.)
Note that this dimensionality reduction tSNE plot does not exactly match the published figure, as the algorithm is stochastic, and a random seed was not set at the time of creating publication figure.
palette_groups_mut <- palette_groups
palette_groups_mut <- c(palette_groups_mut, "EC-NB: MYCN Amp" = "#A575D3")
palette_groups_mut <- c(palette_groups_mut, "EC-NB: FOXR2+, MYCN NonAmp" = "#f382c6")
palette_groups_mut <- c(palette_groups_mut, "EC-NB: FOXR2-, MYCN NonAmp" = "#215ba3")
perplexity <- 30
meta_ecnb <- meta_ecnb %>%
mutate(Group = case_when(
MYCN == "Amp" ~ "EC-NB: MYCN Amp",
MYCN == "NonAmp" & FOXR2_positive == "Y" ~ "EC-NB: FOXR2+, MYCN NonAmp",
MYCN == "NonAmp" & FOXR2_positive == "N" ~ "EC-NB: FOXR2-, MYCN NonAmp"))
meta_all <- bind_rows(meta_bulk, meta_ecnb)
non_stage4_HR_ecnb <- meta_ecnb %>%
filter(Stage != "4" | Risk != "HR" | is.na(Stage) | is.na(Risk)) %>%
.$ID
length(non_stage4_HR_ecnb)
## [1] 387
length(filter(meta_ecnb, !(ID %in% non_stage4_HR_ecnb))$ID)
## [1] 192
set.seed(100)
tsne_calc <- scores_all %>%
select(c("Signature", test_samples, filter(meta_ecnb, !(ID %in% non_stage4_HR_ecnb))$ID)) %>%
tibble::column_to_rownames(var = "Signature") %>%
as.matrix() %>%
ssgsea_dim_red(method = "tsne", return_df = TRUE, perplexity = perplexity)
## @ returning data...
(tsne_plot <- tsne_calc %>% left_join(meta_all, by = c("Sample" = "ID")) %>%
.[which(!(.$Sample %in% non_stage4_HR_ecnb)),] %>%
plot_dim_red(method = "tsne", palette = palette_groups_mut,
point_size = 2, alpha = 0.8, title = glue("High risk & Stage 4 EC-NB only, Perplexity: {perplexity}")) + square_theme + large_text)
Check if distributions are similar:
scores_tidy %>%
ggplot(aes(x = Sample, y = Score)) +
geom_boxplot(aes(fill = Group), outlier.shape = NA) +
scale_fill_manual(values = palette_groups) +
rotate_x() +
ggtitle("Distribution of ssGSEA scores (across signatures) in each sample") +
theme(axis.text.x = element_blank())
Pediatric brain tumor cohort:
scores_tidy %>%
filter(Sample %in% test_samples) %>%
group_by(Sample) %>%
arrange(desc(Score)) %>%
top_n(10, Score) %>%
ggplot(aes(x = Sample, y = Score)) +
geom_point(aes(colour = Class), size = 2, alpha = 0.8) +
scale_colour_manual(values = palette_class) +
facet_grid(~ Group, space = "free", scales = "free_x") +
rotate_x() +
ggtitle("Top 10 signatures in each sample, PBTs, coloured by cell class")
Gartlgruber 2021 extracranial neuroblastoma:
scores_tidy %>%
filter(Sample %in% meta_ecnb$ID) %>%
group_by(Sample) %>%
arrange(desc(Score)) %>%
top_n(10, Score) %>%
ggplot(aes(x = Sample, y = Score)) +
geom_point(aes(colour = Class), size = 0.5, alpha = 0.8) +
scale_colour_manual(values = palette_class) +
facet_grid(~ Group, space = "free", scales = "free_x") +
rotate_x() +
theme(axis.text.x = element_blank()) +
ggtitle("Top 10 signatures in each sample, EC-NB Gartlgruber 2021, coloured by cell class")
Top score in each sample, across pediatric brain tumor groups:
scores_top_pbt <- scores_tidy %>%
filter(Sample %in% test_samples) %>%
group_by(Sample) %>%
arrange(desc(Score)) %>%
top_n(1, Score) %>%
arrange(desc(Score)) %>%
mutate(Sample = factor(Sample, levels = unique(.$Sample)))
scores_top_pbt %>%
ggplot(aes(x = Sample, y = Score)) +
geom_col(aes(fill = Class), width = 0.8) +
scale_fill_manual(values = palette_class) +
xlab("Sample") + ylab(NULL) +
geom_text(aes(label = Signature), angle = 90, nudge_y = 200, hjust = 0) +
facet_grid(~ Group, scales = "free_x", space = "free") +
no_legend() +
theme_min2() +
theme(axis.text.x = element_blank()) +
coord_cartesian(ylim = c(22500, 35000))
# Print top signature for each sample
scores_top_pbt %>% group_by(Group) %>%
mutate(n_group = n()) %>%
group_by(Group, Signature) %>%
summarize(n = n(), prop = n / n_group) %>%
distinct() %>%
arrange(Group, desc(n))
## `summarise()` has grouped output by 'Group', 'Signature'. You can override
## using the `.groups` argument.
Top matching signature per tumor type:
scores_tidy %>%
filter(Sample %in% test_samples) %>%
group_by(Sample) %>%
arrange(desc(Score)) %>%
top_n(1, Score) %>%
mutate(Group = factor(Group,
levels = c("NB-FOXR2",
"ETMR",
"MB-SHH",
"MB-WNT",
"DIPG-H3K27M-FOXR2",
"DIPG-H3K27M",
"HGG-IDH",
"EP-PFA",
"HGG-H3.3G34R/V"))) %>%
ggplot(aes(x = Group)) +
geom_bar(aes(fill = Class), size = 2, alpha = 0.8) +
scale_fill_manual(values = palette_class) +
theme(axis.text.x = element_text(angle = 30, vjust = 1, hjust=1))
Top signature in NB-FOXR2:
scores_top_nbfoxr2 <- scores_top_pbt %>% filter(Group == "NB-FOXR2")
scores_top_nbfoxr2 %>%
ungroup() %>%
count(Signature, Class) %>%
ggplot(aes(x = reorder(Signature, n),
y = n,
fill = Class)) +
geom_col() +
scale_fill_manual(values = palette_class) +
coord_flip() +
xlab("Signature") +
ylab("Number of samples") +
large_text +
theme(legend.position = "bottom")
Top signature in extracranial neuroblastoma (Gartlgruber cohort):
scores_top_ec <- scores_tidy %>%
filter(Sample %in% meta_ecnb$ID) %>%
group_by(Sample) %>%
arrange(desc(Score)) %>%
top_n(1, Score) %>%
arrange(desc(Score)) %>%
mutate(Sample = factor(Sample, levels = unique(.$Sample)))
scores_top_ec %>%
ggplot(aes(x = Sample, y = Score)) +
geom_col(aes(fill = Class, color = Class)) +
scale_color_manual(values = palette_class) +
scale_fill_manual(values = palette_class) +
xlab("Sample") + ylab(NULL) +
theme(axis.text.x = element_blank()) +
coord_cartesian(ylim = c(22500, 30000))
# Print percentage of EC-NB samples mapping to each signature
(prop <- scores_top_ec %>%
group_by(Signature) %>%
count() %>%
mutate(prop = n/sum(.$n)) %>%
arrange(desc(n)))
# Print total number of EC-NB samples in Gartlgruber cohort
sum(prop$n)
## [1] 579
Save top scoring signatures for each tumor in a supplementary table:
TABLE_top_scores <- bind_rows(scores_top_pbt, scores_top_ec) %>%
select(Sample, Group, Signature, Score, Reference_sample,
Reference_age, Cell_type, Dataset, Species, Class) %>%
arrange(Group)
dim(TABLE_top_scores)
## [1] 700 10
table(as.character(TABLE_top_scores$Group))
##
## DIPG-H3K27M DIPG-H3K27M-FOXR2 EC-NB EP-PFA
## 19 6 579 13
## ETMR HGG-H3.3G34R/V HGG-IDH MB-SHH
## 12 18 10 8
## MB-WNT NB-FOXR2
## 10 25
rr_write_tsv(TABLE_top_scores,
glue("{out}/TABLE_top_scores.tsv"),
"Top scoring ssGSEA signature for each bulk tumor")
## ...writing description of TABLE_top_scores.tsv to public/output/06.1/TABLE_top_scores.desc
Define a helper function for per-sample plots of the top 10 signatures:
# Extract top 10 per sample
top_scores <- scores_tidy %>%
group_by(Sample) %>%
arrange(desc(Score)) %>%
top_n(10, Score)
#' @param goi Character, group of interest
plot_scores_per_sample <- function(goi) {
samples <- top_scores %>% filter(Group == goi) %>% pull(Sample) %>% unique()
p1 <- map(samples,
~ top_scores %>%
filter(Sample == .x) %>%
arrange(Score) %>%
mutate(Signature = factor(Signature, levels = .$Signature)) %>%
ggplot(aes(x = Signature, y = Score)) +
geom_col(aes(fill = Class), width = 0.5) +
scale_fill_manual(values = palette_class) +
theme_min2() +
rotate_x() +
xlab(NULL) + ylab(NULL) +
ggtitle(.x) +
no_legend() +
coord_flip(ylim = c(20000, 30000)))
plot_grid(plotlist = p1, align = "hv", axis = "rltb", ncol = 4)
}
Plot top 10 signatures per NB-FOXR2 sample:
plot_scores_per_sample("NB-FOXR2")
In this section, we perform a comparative analysis of ssGSEA enrichment scores between tumor types. We'll broadly divide brain tumors into neuronal and glial tumor types, and compare NB-FOXR2 to each of these broad groups. Here, we exclude H3.3G34R/V tumors, since this tumor has an interneuron signal (Chen et al, Cell 2020).
# Label tumors as neuronal or glial, remove G34
scores_tidy_comp <- scores_tidy %>%
filter(!grepl("EC", Group)) %>%
mutate(Group_comp = case_when(
Group == "NB-FOXR2" ~ "NB-FOXR2",
grepl("MB|ETMR|HGNET", Group) ~ "Neuronal",
grepl("HGG|PFA|DIPG", Group) &
Group != "HGG-H3.3G34R/V" ~ "Glial",
))
# sanity check
any(scores_tidy$Signature %in% sigs_indiscriminant)
## [1] FALSE
# compute the mean per signature per group
scores_group_means <- scores_tidy_comp %>%
mutate(Score = scales::rescale(Score, to = c(0, 1))) %>%
group_by(Signature, Group_comp) %>%
summarise(group_mean = mean(Score))
## `summarise()` has grouped output by 'Signature'. You can override using the
## `.groups` argument.
# using broom to run pairwise t tests between groups based on the score
# NOTE: see below on differences in defaults between pairwise.t.test and t.test
# https://stackoverflow.com/a/11457871
scores_tidy_ttest <- scores_tidy_comp %>%
group_by(Signature) %>%
summarise(ttest = list(pairwise.t.test(Score, Group_comp, p.adjust = "none",
paired = FALSE,
pool.sd = FALSE))) %>%
mutate(ttest = map(ttest, broom::tidy)) %>%
unnest(cols = c(ttest)) %>%
mutate(p_adj = p.adjust(p.value, method = "BH"))
# add the group estimates after pairwise tests for later analysis
scores_tidy_ttest <- scores_tidy_ttest %>%
left_join(scores_group_means, by = c("Signature", "group1" = "Group_comp")) %>%
dplyr::rename(group1_mean = group_mean) %>%
left_join(scores_group_means, by = c("Signature", "group2" = "Group_comp")) %>%
dplyr::rename(group2_mean = group_mean)
write_tsv(scores_tidy_ttest, glue("{out}/ssgsea_scores_tests.tsv"))
Data wrangling to make sure the fold-changes are in the same direction for each comparison:
scores_tidy_ttest_foxr2 <- scores_tidy_ttest %>%
filter(group1 == "NB-FOXR2" | group2 == "NB-FOXR2") %>%
mutate(log10p = -log10(p_adj),
log2fc = log2(group1_mean / group2_mean)) %>%
# put all group1 as NB-FOXR2, and reverse L2FC if needed
# this is because the pairwise t-tests do not perform two tests for each
# pair of groups (i.e. only A vs B, not A vs B and B vs A)
mutate(
group_tmp = group1,
# if the Group2 is NB-FOXR2, then the LFC sign needs to be switched
log2fc = ifelse(group2 == "NB-FOXR2", -log2fc, log2fc),
# if the Group2 is NB-FOXR2, swap group1 and 2
group1 = ifelse(group2 == "NB-FOXR2", "NB-FOXR2", group1),
group2 = ifelse(group2 == "NB-FOXR2", group_tmp, group2)
) %>%
select(-group_tmp) %>%
mutate(signif = ifelse(p_adj < 0.01, TRUE, FALSE)) %>%
left_join(cell_type_anno_all, by = c("Signature" = "Cluster"))
Visualize the results as volcano plots:
df1 <- scores_tidy_ttest_foxr2 %>%
filter(group2 == "Neuronal")
p1 <- df1 %>%
mutate(Class = ifelse(log2fc < 0 | p_adj > 0.01, "Other", Class)) %>%
ggplot(aes(x = log2fc, y = -log10(p_adj))) +
geom_hline(yintercept = 2, color = "gray50", linetype = 2) +
geom_vline(xintercept = 0, color = "gray50", linetype = 2) +
geom_point(aes(color = Class, size = -log10(p_adj)), alpha = 0.8) +
scale_color_manual(values = palette_class) +
no_legend() +
ggtitle("NB-FOXR2 vs neuronal tumors") +
coord_cartesian(xlim = c(-0.3, 0.3))
df2 <- scores_tidy_ttest_foxr2 %>%
filter(group2 == "Glial")
p2 <- df2 %>%
mutate(Class = ifelse(log2fc < 0 | p_adj > 0.01, "Other", Class)) %>%
ggplot(aes(x = log2fc, y = -log10(p_adj))) +
geom_hline(yintercept = 2, color = "gray50", linetype = 2) +
geom_vline(xintercept = 0, color = "gray50", linetype = 2) +
geom_point(aes(color = Class, size = -log10(p_adj)), alpha = 0.8) +
scale_color_manual(values = palette_class) +
no_legend() +
ggtitle("NB-FOXR2 vs glial tumors") +
coord_cartesian(xlim = c(-0.3, 0.3))
plot_grid(p1, p2, nrow = 1)
# with signature labels
p1 <- p1 + geom_text_repel(data = df1 %>% filter(-log10(p_adj) > 5 & log2fc > 0),
aes(label = Signature, color = Class),
size = 4)
p2 <- p2 + geom_text_repel(data = df2 %>% filter(-log10(p_adj) > 10 & log2fc > 0),
aes(label = Signature, color = Class),
size = 4)
plot_grid(p1, p2, nrow = 1)
This document was last rendered on:
## 2024-11-05 10:36:03
The git repository and last commit:
## Local: main /project/kleinman/bhavyaa.chandarana/from_hydra/2023-05-NB-FOXR2/public
## Remote: main @ origin (https://github.com/fungenomics/NB-FOXR2.git)
## Head: [72d1c5c] 2024-11-04: Add DOI badge
The random seed was set with set.seed(100)
The R session info:
## ─ Session info ───────────────────────────────────────────────────────────────
## setting value
## version R version 4.1.2 (2021-11-01)
## os Rocky Linux 8.10 (Green Obsidian)
## system x86_64, linux-gnu
## ui X11
## language (EN)
## collate en_US.UTF-8
## ctype en_US.UTF-8
## tz America/Toronto
## date 2024-11-05
## pandoc 1.19.2.1 @ /cvmfs/soft.computecanada.ca/gentoo/2020/usr/bin/ (via rmarkdown)
##
## ─ Packages ───────────────────────────────────────────────────────────────────
## ! package * version date (UTC) lib source
## P abind 1.4-5 2016-07-21 [?] CRAN (R 4.1.2)
## P backports 1.4.1 2021-12-13 [?] CRAN (R 4.1.2)
## P BiocGenerics 0.40.0 2021-10-26 [?] Bioconductor
## P BiocManager 1.30.15 2021-05-11 [?] CRAN (R 4.1.2)
## P bit 4.0.4 2020-08-04 [?] CRAN (R 4.1.2)
## P bit64 4.0.5 2020-08-30 [?] CRAN (R 4.1.2)
## P broom * 1.0.1 2022-08-29 [?] CRAN (R 4.1.2)
## P bslib 0.3.1 2021-10-06 [?] CRAN (R 4.1.2)
## P cachem 1.0.6 2021-08-19 [?] CRAN (R 4.1.2)
## P callr 3.7.6 2024-03-25 [?] RSPM
## P cellranger 1.1.0 2016-07-27 [?] CRAN (R 4.1.2)
## P circlize 0.4.15 2022-05-10 [?] CRAN (R 4.1.2)
## P cli 3.6.1 2023-03-23 [?] RSPM (R 4.1.2)
## P clue 0.3-64 2023-01-31 [?] CRAN (R 4.1.2)
## P cluster 2.1.2 2021-04-17 [?] CRAN (R 4.1.2)
## P codetools 0.2-18 2020-11-04 [?] CRAN (R 4.1.2)
## P colorspace 2.0-2 2021-06-24 [?] CRAN (R 4.1.2)
## P ComplexHeatmap * 2.10.0 2021-10-26 [?] Bioconductor
## P conflicted 1.2.0 2023-02-01 [?] CRAN (R 4.1.2)
## P cowplot * 1.1.1 2020-12-30 [?] CRAN (R 4.1.2)
## P crayon 1.4.2 2021-10-29 [?] CRAN (R 4.1.2)
## P crosstalk 1.2.0 2021-11-04 [?] CRAN (R 4.1.2)
## P data.table * 1.14.2 2021-09-27 [?] CRAN (R 4.1.2)
## P deldir 1.0-6 2021-10-23 [?] CRAN (R 4.1.2)
## P devtools 2.4.5 2022-10-11 [?] CRAN (R 4.1.2)
## P digest 0.6.35 2024-03-11 [?] CRAN (R 4.1.2)
## P doParallel 1.0.16 2020-10-16 [?] CRAN (R 4.1.2)
## P dplyr * 1.1.1 2023-03-22 [?] CRAN (R 4.1.2)
## P DT 0.30 2023-10-05 [?] CRAN (R 4.1.2)
## P ellipsis 0.3.2 2021-04-29 [?] CRAN (R 4.1.2)
## P evaluate 0.23 2023-11-01 [?] CRAN (R 4.1.2)
## P fansi 1.0.2 2022-01-14 [?] CRAN (R 4.1.2)
## P farver 2.1.0 2021-02-28 [?] CRAN (R 4.1.2)
## P fastmap 1.1.0 2021-01-25 [?] CRAN (R 4.1.2)
## P fitdistrplus 1.1-6 2021-09-28 [?] CRAN (R 4.1.2)
## P foreach 1.5.1 2020-10-15 [?] CRAN (R 4.1.2)
## P fs 1.5.2 2021-12-08 [?] CRAN (R 4.1.2)
## P future 1.25.0 2022-04-24 [?] CRAN (R 4.1.2)
## P future.apply 1.8.1 2021-08-10 [?] CRAN (R 4.1.2)
## P generics 0.1.3 2022-07-05 [?] CRAN (R 4.1.2)
## P GetoptLong 1.0.5 2020-12-15 [?] CRAN (R 4.1.2)
## P ggExtra * 0.10.0 2022-03-23 [?] CRAN (R 4.1.2)
## P ggplot2 * 3.4.2 2023-04-03 [?] CRAN (R 4.1.2)
## P ggrepel * 0.9.1 2021-01-15 [?] CRAN (R 4.1.2)
## P ggridges 0.5.3 2021-01-08 [?] CRAN (R 4.1.2)
## P git2r 0.29.0 2021-11-22 [?] CRAN (R 4.1.2)
## P GlobalOptions 0.1.2 2020-06-10 [?] CRAN (R 4.1.2)
## P globals 0.14.0 2020-11-22 [?] CRAN (R 4.1.2)
## P glue * 1.6.2 2022-02-24 [?] CRAN (R 4.1.2)
## P goftest 1.2-3 2021-10-07 [?] CRAN (R 4.1.2)
## P gridExtra 2.3 2017-09-09 [?] CRAN (R 4.1.2)
## P gtable 0.3.0 2019-03-25 [?] CRAN (R 4.1.2)
## P here * 1.0.1 2020-12-13 [?] CRAN (R 4.1.2)
## P highr 0.9 2021-04-16 [?] CRAN (R 4.1.2)
## P hms 1.1.1 2021-09-26 [?] CRAN (R 4.1.2)
## P htmltools 0.5.2 2021-08-25 [?] CRAN (R 4.1.2)
## P htmlwidgets 1.5.4 2021-09-08 [?] CRAN (R 4.1.2)
## P httpuv 1.6.5 2022-01-05 [?] CRAN (R 4.1.2)
## P httr 1.4.2 2020-07-20 [?] CRAN (R 4.1.2)
## P ica 1.0-2 2018-05-24 [?] CRAN (R 4.1.2)
## P igraph 2.0.3 2024-03-13 [?] CRAN (R 4.1.2)
## P IRanges 2.28.0 2021-10-26 [?] Bioconductor
## P irlba 2.3.5 2021-12-06 [?] CRAN (R 4.1.2)
## P iterators 1.0.13 2020-10-15 [?] CRAN (R 4.1.2)
## P jquerylib 0.1.4 2021-04-26 [?] CRAN (R 4.1.2)
## P jsonlite 1.8.8 2023-12-04 [?] CRAN (R 4.1.2)
## P KernSmooth 2.23-20 2021-05-03 [?] CRAN (R 4.1.2)
## P knitr 1.37 2021-12-16 [?] CRAN (R 4.1.2)
## P labeling 0.4.2 2020-10-20 [?] CRAN (R 4.1.2)
## P later 1.3.0 2021-08-18 [?] CRAN (R 4.1.2)
## P lattice 0.20-45 2021-09-22 [?] CRAN (R 4.1.2)
## P lazyeval 0.2.2 2019-03-15 [?] CRAN (R 4.1.2)
## P leiden 0.3.9 2021-07-27 [?] CRAN (R 4.1.2)
## P lifecycle 1.0.3 2022-10-07 [?] CRAN (R 4.1.2)
## P listenv 0.8.0 2019-12-05 [?] CRAN (R 4.1.2)
## P lmtest 0.9-39 2021-11-07 [?] CRAN (R 4.1.2)
## P magrittr * 2.0.3 2022-03-30 [?] CRAN (R 4.1.2)
## P MASS 7.3-54 2021-05-03 [?] CRAN (R 4.1.2)
## P Matrix 1.3-4 2021-06-01 [?] CRAN (R 4.1.2)
## P matrixStats 0.61.0 2021-09-17 [?] CRAN (R 4.1.2)
## P memoise 2.0.1 2021-11-26 [?] CRAN (R 4.1.2)
## P mgcv 1.8-38 2021-10-06 [?] CRAN (R 4.1.2)
## P mime 0.12 2021-09-28 [?] CRAN (R 4.1.2)
## P miniUI 0.1.1.1 2018-05-18 [?] CRAN (R 4.1.2)
## P munsell 0.5.0 2018-06-12 [?] CRAN (R 4.1.2)
## P nlme 3.1-153 2021-09-07 [?] CRAN (R 4.1.2)
## P parallelly 1.30.0 2021-12-17 [?] CRAN (R 4.1.2)
## P patchwork 1.1.1 2020-12-17 [?] CRAN (R 4.1.2)
## P pbapply * 1.5-0 2021-09-16 [?] CRAN (R 4.1.2)
## P pillar 1.9.0 2023-03-22 [?] RSPM (R 4.1.2)
## P pkgbuild 1.4.2 2023-06-26 [?] CRAN (R 4.1.2)
## P pkgconfig 2.0.3 2019-09-22 [?] CRAN (R 4.1.2)
## P pkgload 1.3.3 2023-09-22 [?] CRAN (R 4.1.2)
## P plotly 4.10.0 2021-10-09 [?] CRAN (R 4.1.2)
## P plyr 1.8.6 2020-03-03 [?] CRAN (R 4.1.2)
## P png 0.1-7 2013-12-03 [?] CRAN (R 4.1.2)
## P polyclip 1.10-0 2019-03-14 [?] CRAN (R 4.1.2)
## P prettyunits 1.1.1 2020-01-24 [?] CRAN (R 4.1.2)
## P processx 3.8.4 2024-03-16 [?] RSPM
## P profvis 0.3.8 2023-05-02 [?] CRAN (R 4.1.2)
## P promises 1.2.0.1 2021-02-11 [?] CRAN (R 4.1.2)
## P ps 1.7.6 2024-01-18 [?] RSPM
## P purrr * 1.0.1 2023-01-10 [?] CRAN (R 4.1.2)
## P R6 2.5.1 2021-08-19 [?] CRAN (R 4.1.2)
## P RANN 2.6.1 2019-01-08 [?] CRAN (R 4.1.2)
## P RColorBrewer * 1.1-2 2014-12-07 [?] CRAN (R 4.1.2)
## P Rcpp 1.0.8 2022-01-13 [?] CRAN (R 4.1.2)
## P RcppAnnoy 0.0.19 2021-07-30 [?] CRAN (R 4.1.2)
## P readr * 2.1.1 2021-11-30 [?] CRAN (R 4.1.2)
## P readxl * 1.3.1 2019-03-13 [?] CRAN (R 4.1.2)
## P remotes 2.4.2.1 2023-07-18 [?] CRAN (R 4.1.2)
## P renv 1.0.3 2023-09-19 [?] CRAN (R 4.1.2)
## P reshape2 1.4.4 2020-04-09 [?] CRAN (R 4.1.2)
## P reticulate 1.23 2022-01-14 [?] CRAN (R 4.1.2)
## P rjson 0.2.21 2022-01-09 [?] CRAN (R 4.1.2)
## P rlang 1.1.3 2024-01-10 [?] CRAN (R 4.1.2)
## P rmarkdown 2.11 2021-09-14 [?] CRAN (R 4.1.2)
## P ROCR 1.0-11 2020-05-02 [?] CRAN (R 4.1.2)
## P rpart 4.1-15 2019-04-12 [?] CRAN (R 4.1.2)
## P rprojroot 2.0.2 2020-11-15 [?] CRAN (R 4.1.2)
## P Rtsne * 0.15 2018-11-10 [?] CRAN (R 4.1.2)
## P S4Vectors 0.32.4 2022-03-24 [?] Bioconductor
## P sass 0.4.0 2021-05-12 [?] CRAN (R 4.1.2)
## P scales 1.2.1 2022-08-20 [?] CRAN (R 4.1.2)
## P scattermore 0.7 2020-11-24 [?] CRAN (R 4.1.2)
## P sctransform 0.3.3 2022-01-13 [?] CRAN (R 4.1.2)
## P sessioninfo 1.2.2 2021-12-06 [?] CRAN (R 4.1.2)
## P Seurat * 4.0.0 2021-01-30 [?] CRAN (R 4.1.2)
## P SeuratObject * 4.0.4 2021-11-23 [?] CRAN (R 4.1.2)
## P shape 1.4.6 2021-05-19 [?] CRAN (R 4.1.2)
## P shiny 1.7.1 2021-10-02 [?] CRAN (R 4.1.2)
## P spatstat 1.64-1 2020-05-12 [?] CRAN (R 4.1.2)
## P spatstat.data 2.1-2 2021-12-17 [?] CRAN (R 4.1.2)
## P spatstat.utils 2.3-0 2021-12-12 [?] CRAN (R 4.1.2)
## P stringi 1.7.6 2021-11-29 [?] CRAN (R 4.1.2)
## P stringr * 1.5.0 2022-12-02 [?] CRAN (R 4.1.2)
## P survival 3.2-13 2021-08-24 [?] CRAN (R 4.1.2)
## P tensor 1.5 2012-05-05 [?] CRAN (R 4.1.2)
## P tibble * 3.2.1 2023-03-20 [?] RSPM (R 4.1.2)
## P tidyr * 1.3.0 2023-01-24 [?] CRAN (R 4.1.2)
## P tidyselect 1.2.0 2022-10-10 [?] CRAN (R 4.1.2)
## P tzdb 0.3.0 2022-03-28 [?] CRAN (R 4.1.2)
## P urlchecker 1.0.1 2021-11-30 [?] CRAN (R 4.1.2)
## P usethis 2.2.2 2023-07-06 [?] CRAN (R 4.1.2)
## P utf8 1.2.2 2021-07-24 [?] CRAN (R 4.1.2)
## P uwot 0.1.11 2021-12-02 [?] CRAN (R 4.1.2)
## P vctrs 0.6.5 2023-12-01 [?] CRAN (R 4.1.2)
## P viridis * 0.5.1 2018-03-29 [?] RSPM (R 4.1.2)
## P viridisLite * 0.3.0 2018-02-01 [?] CRAN (R 4.1.2)
## P vroom 1.5.7 2021-11-30 [?] CRAN (R 4.1.2)
## P withr 2.5.0 2022-03-03 [?] CRAN (R 4.1.2)
## P xfun 0.29 2021-12-14 [?] CRAN (R 4.1.2)
## P xtable 1.8-4 2019-04-21 [?] CRAN (R 4.1.2)
## P yaml 2.2.1 2020-02-01 [?] CRAN (R 4.1.2)
## P zoo 1.8-9 2021-03-09 [?] CRAN (R 4.1.2)
##
## [1] /project/kleinman/bhavyaa.chandarana/from_hydra/2023-05-NB-FOXR2/public/renv/library/R-4.1/x86_64-pc-linux-gnu
## [2] /home/kleinman/bhavyaa.chandarana/.cache/R/renv/sandbox/R-4.1/x86_64-pc-linux-gnu/145cef2c
##
## P ── Loaded and on-disk path mismatch.
##
## ──────────────────────────────────────────────────────────────────────────────
A project of the Kleinman Lab at McGill University, using the rr reproducible research template.