Splitting up #921: Immune-deconv changes

analyses/immune-deconv/01-immune-deconv.R

@@ -3,6 +3,9 @@
 # Function:
 # Script to perform immune characterization using xCell etc.
 
+# Find the root directory of this repository
+root_dir <- rprojroot::find_root(rprojroot::has_dir(".git"))
+
 # load libraries
 suppressPackageStartupMessages(library(optparse))
 suppressPackageStartupMessages(library(tidyverse))
@@ -73,10 +76,19 @@ deconv <- function(expr.input, method) {
   # get data
   expr.input <- get(expr.input)
 
+  # Import standard color palettes for project
+  histology_label_mapping <- readr::read_tsv(
+    file.path(root_dir, "figures", "palettes", "histology_label_color_table.tsv")) %>%
+    # Select just the columns we will need for plotting
+    dplyr::select(Kids_First_Biospecimen_ID, display_group, display_order, hex_codes) %>%
+    # Reorder display_group based on display_order
+    dplyr::mutate(display_group = forcats::fct_reorder(display_group, display_order))
+
   # subset clinical
   clin.sub  <- clin %>%
     filter(Kids_First_Biospecimen_ID %in% colnames(expr.input)) %>%
-    dplyr::select(Kids_First_Biospecimen_ID, broad_histology, short_histology, molecular_subtype)
+    dplyr::inner_join(histology_label_mapping, by = "Kids_First_Biospecimen_ID") %>%
+    dplyr::select(Kids_First_Biospecimen_ID, broad_histology, display_group, molecular_subtype)
 
   # deconvolute using specified method
   res <- deconvolute(gene_expression = as.matrix(expr.input), method = method)

analyses/immune-deconv/02-summary-plots.R

@@ -14,7 +14,7 @@ source(file.path(root_dir, "analyses", "immune-deconv",
                  "util", "pubTheme.R"))
 
 option_list <- list(
-  make_option(c("-i", "--input"), type = "character", help = "Immunedeconv output from 01-immune.deconv.R (.RData)"), 
+  make_option(c("-i", "--input"), type = "character", help = "Immunedeconv output from 01-immune.deconv.R (.RData)"),
   make_option(c("-o",  "--output"), type = "character", help = "Output directory")
 )
 
@@ -29,13 +29,9 @@ deconvout <- opt$input
 output <- opt$output
 load(deconvout) 
 
-# if short histology is Medulloblastoma or ATRT, use them as broad histology
-deconv.res$broad_histology <- ifelse(deconv.res$short_histology %in% c("ATRT", "Medulloblastoma"),
-                                     deconv.res$short_histology, deconv.res$broad_histology)
-
 # add molecular subtype info
-deconv.res$broad_histology <- ifelse(is.na(deconv.res$molecular_subtype), 
-                                     deconv.res$broad_histology, 
+deconv.res$broad_histology <- ifelse(is.na(deconv.res$molecular_subtype),
+                                     deconv.res$broad_histology,
                                      paste0(deconv.res$broad_histology, '-', deconv.res$molecular_subtype))
 
 
@@ -53,35 +49,35 @@ method2 <- deconv.res %>%
 # first, define a function to create heatmap of
 # average immune scores per histology per cell type
 create.heatmap <- function(deconv.method, title, fileout) {
-  
+
   # assign labels
   non.brain.tumors <- c("Histiocytic tumor", "Lymphomas")
-  
+
   # create labels: count of samples per histology
   annot <- deconv.method %>%
     dplyr::select(broad_histology, sample) %>%
     unique() %>%
     group_by(broad_histology) %>%
     summarise(label = n()) %>%
     mutate(label = paste0(broad_histology,' (',label,')'))
-  
+
   # add labels to actual data
   deconv.method <- merge(deconv.method, annot, by = 'broad_histology')
-  
+
   # calculate average scores per cell type per histology
-  deconv.method <- deconv.method %>% 
+  deconv.method <- deconv.method %>%
     filter(!cell_type %in% c("microenvironment score", "stroma score", "immune score")) %>%
     group_by(cell_type, label) %>%
     dplyr::summarise(mean = mean(fraction)) %>%
     # convert into matrix of cell type vs histology
-    spread(key = label, value = mean) %>% 
+    spread(key = label, value = mean) %>%
     column_to_rownames('cell_type')
-  
+
   # plot non-brain and brain tumors separately
-  pdf(file = fileout, width = 15, height = 10)
+  pdf(file = fileout, width = 15, height = 15)
   # non-brain tumors
-  mat <- deconv.method %>% 
-    dplyr::select(grep(paste0(non.brain.tumors, collapse="|"), colnames(deconv.method), value = TRUE)) 
+  mat <- deconv.method %>%
+    dplyr::select(grep(paste0(non.brain.tumors, collapse="|"), colnames(deconv.method), value = TRUE))
   if(ncol(mat) > 1){
     mat <- mat %>%
       rownames_to_column('celltype') %>%
@@ -90,11 +86,11 @@ create.heatmap <- function(deconv.method, title, fileout) {
       t() %>%
       pheatmap(fontsize = 10,
                scale = "column", angle_col = 45,
-               main = "Average immune scores normalized by rows\nNon-Brain Tumors", 
+               main = "Average immune scores normalized by rows\nNon-Brain Tumors",
                annotation_legend = T, cellwidth = 15, cellheight = 15)
   }
-  
-  # brain tumors 
+
+  # brain tumors
   mat <- deconv.method %>%
     dplyr::select(grep(paste0(non.brain.tumors, collapse="|"), colnames(deconv.method), invert = TRUE, value = TRUE))
   if(ncol(mat) > 1){
@@ -103,18 +99,18 @@ create.heatmap <- function(deconv.method, title, fileout) {
       filter_at(vars(-celltype), any_vars(. != 0)) %>%
       column_to_rownames('celltype') %>%
       t() %>%
-      pheatmap(fontsize = 10, 
+      pheatmap(fontsize = 10,
                scale = "column", angle_col = 45,
-               main = "Average immune scores normalized by rows\nBrain Tumors", 
+               main = "Average immune scores normalized by rows\nBrain Tumors",
                annotation_legend = T, cellwidth = 15, cellheight = 15)
   }
-  
+
   dev.off()
 }
 
 # next, plot a correlation heatmap between xCell and the second specified method
 # only take common cell types between both methods
-common.types <- intersect(method1$cell_type, method2$cell_type) 
+common.types <- intersect(method1$cell_type, method2$cell_type)
 method1.sub <- method1 %>%
   filter(cell_type %in% common.types) %>%
   mutate(!!method1.name := fraction) %>%
@@ -141,22 +137,22 @@ total.labels <- total %>%
 total <- merge(total, total.labels, by = 'broad_histology')
 
 # calculate correlation per cell type per histology
-total <- total %>% 
+total <- total %>%
   group_by(cell_type, label) %>%
   dplyr::summarise(corr = cor(!!sym(method1.name), !!sym(method2.name))) %>%
-  spread(key = label, value = corr) %>% 
+  spread(key = label, value = corr) %>%
   column_to_rownames('cell_type') %>%
   replace(is.na(.), 0)
 
 # replace space from method names for output filename
-m1 <- gsub(" ","",method1.name) 
+m1 <- gsub(" ","",method1.name)
 m2 <- gsub(" ","",method2.name)
 
 # create correlation plot for overlapping cell types between both methods
-pdf(file = file.path(output, paste0("corrplot_", m1, "_vs_", m2, ".pdf")), 
-    width = 16, height = 8)
+pdf(file = file.path(output, paste0("corrplot_", m1, "_vs_", m2, ".pdf")),
+    width = 16, height = 20)
 corrplot(t(total), method = "circle", type = 'full', win.asp = 0.5,
-         addCoef.col = "#888888", number.cex = .7,  
+         addCoef.col = "#888888", number.cex = .7,
          tl.col = "black", number.font = 2,
          is.corr = FALSE, tl.cex = 0.8,
          mar = c(0, 0, 0, 5),