refactor: trajectory_pseudotime_heatmap change back to script

scripts/plotting/trajectory_pseudotime_heatmap.R

@@ -1,4 +1,3 @@
-trajectory_pseudotime_heatmap <- function(sce, clustering, lineage_num=1, n_top_genes=50, out_file=NULL) {
 
 library(Seurat)
 library(slingshot)
@@ -10,23 +9,36 @@ library(tibble)
 library(pheatmap)
 library(SingleCellExperiment)
 
-# Pre-load some inputs
-# sce <- readRDS('/Users/jonrob/Documents/NBIS/LTS_projects/d_angeletti_1910/analysis/trajectory_slingshot/trajectory_07/sce_object.rds')
-# DATA <- readRDS('/Users/jonrob/Documents/NBIS/LTS_projects/d_angeletti_1910/analysis/06_cluster/seurat_object.rds')
-# DATA <- subset(DATA, cells=sce$slingshot@rownames)
-# clustering <- as.factor(DATA@meta.data[, 'HC_16']) %>% setNames(colnames(DATA))
+################################
+### DEFINE PARAMETERS & DATA ###
+################################
+# specify params
+lineage_num <- 1
+n_top_genes <- 50
+
+# specify directories
+project_dir <- '/Users/jonrob/Documents/NBIS/LTS_projects/d_angeletti_1910/scMouseBcellFlu'
+trajectory_dir <- paste0(project_dir, '/analysis/trajectory_slingshot/trajectory_07')
+
+# specify SCE file (from trajectory analysis) and output filename to which heatmap will be saved
+sce <- readRDS(paste0(trajectory_dir, '/sce_object.rds'))
+out_file <- paste0(trajectory_dir, '/pseudotime_heatmap_curve1.svg')
+
+# define clustering
+DATA <- readRDS(paste0(project_dir, '/analysis/06_cluster/seurat_object.rds'))
+DATA <- subset(DATA, cells=sce$slingshot@rownames)
+clustering <- as.factor(DATA@meta.data[, 'HC_16']) %>% setNames(colnames(DATA))
+rm(DATA)
+#------------------------------
+
 
-
 # function for generating color hues (mimics default ggplot palette)
 gg_color_hue <- function(n) {
   hues = seq(15, 375, length = n + 1)
   hcl(h = hues, l = 65, c = 100)[1:n]
 }
 
 # pre-process clustering data
-# if (is.factor(clustering)) {
-#   clustering <- as.numeric(as.character(clustering)) %>% setNames(names(clustering))
-# }
 clustering <- data.frame(cluster=clustering) %>% rownames_to_column(var='cell')
 if (!any(clustering$cell %in% sce$slingshot@rownames)) {
   stop('Names of clustering input should be cell IDs found in the SCE input.')
@@ -82,14 +94,7 @@ curves_interp <- apply(loess_curves[,-1], 2, function(x) approx(loess_curves$pse
 colnames(curves_interp) <- seq(ncol(curves_interp))
 
 
-
-# # generate heatmap of gene expression as a function of pseudotime
-# heatmap(curves_interp, Colv=NA, distfun=function(x) as.dist(1-cor(t(x))),
-#         labCol=NA, labRow=rownames(curves_interp), col=magma(100,direction=-1),
-#         keep.dendro=F)
-
-
-# terribly inefficient way to create vector to label x-axis
+# inefficient but effective way to create vector to label x-axis
 x_ticks <- seq(0, 1, 0.1)
 keep <- unlist(lapply(x_ticks , function(i) which.min(abs(seq(min(x_ticks), max(x_ticks), length.out=ncol(curves_interp)) - i))))
 xvals <- rep("", ncol(curves_interp))
@@ -98,10 +103,7 @@ xvals[keep] <- x_ticks
 # annotation colors for clusters
 annColors <- list(cluster = gg_color_hue(nlevels(clusters_interp$cluster)) %>% setNames(levels(clusters_interp$cluster)))
 
-# curves_scaled <- t(scale(t(curves_interp), center=F, scale=T))
-
-# png('test_heatmap.png', units='mm', width=200, height=200, res=150)
-# svg('test_heatmap.svg', width=7, height=7)
+# generate heatmap
 p <- pheatmap(curves_interp,
               cluster_cols=F,
               scale='row',
@@ -116,6 +118,7 @@ p <- pheatmap(curves_interp,
               annotation_colors=annColors,
               annotation_names_col=F)
 
+# export heatmap to file
 if (!is.null(out_file)) {
   if (endsWith(casefold(out_file), '.svg')) {
     svg(out_file, width=7, height=7)
@@ -128,12 +131,6 @@ if (!is.null(out_file)) {
   }
 }
 
-return(p)
-
-
-}
-
-
 
 
 
@@ -149,5 +146,3 @@ return(p)
 #          ColSideColors=as.character(clusters_interp)
 #          )
 
-
-