Minimum example to run Monocle2 with either the ICA or a custom function

hsmm.R

# load libs
library(monocle)
library(HSMMSingleCell)
library(rsvd)
library(ggplot2)
library(ggpubr)

# load data
data("HSMM_expr_matrix")
data("HSMM_sample_sheet")
data("HSMM_gene_annotation")

# 1/5 prepare Monocle2 object
pd <- new("AnnotatedDataFrame", data = HSMM_sample_sheet)
fd <- new("AnnotatedDataFrame", data = HSMM_gene_annotation)
HSMM <- newCellDataSet(as.matrix(HSMM_expr_matrix),
                       phenoData = pd, featureData = fd)
HSMM <- estimateSizeFactors(HSMM)
HSMM <- estimateDispersions(HSMM)

# 2/5 select genes for ordering cells
HSMM <- detectGenes(HSMM, min_expr = 0.1)
disp_table <- dispersionTable(HSMM)
ordering_genes <- subset(
  disp_table,
  mean_expression >= 0.5 &
    dispersion_empirical >= 1 * dispersion_fit)$gene_id

HSMM <- setOrderingFilter(HSMM, ordering_genes)
table(fData(HSMM)$use_for_ordering)

# 3/5 run dimension reduction
# option-1: using the default ICA
HSMM_ica <- reduceDimension(
  HSMM,
  norm_method = 'log',
  max_components = 2,
  reduction_method = 'ICA',
  verbose = FALSE)

# option-2: using a custom method, here randomized PCA
randomized_pca_monocle2 <- function(mat){
  # @param tmat A scaled (de-mean & de-var) matrix with features by samples
  # @return A matrix with features by samples projected on PCA space
  rpca_obj <- rpca(t(mat), k=2, center=F, scale=F, retx=T)
  t(rpca_obj$x)
}

HSMM_rpca <- reduceDimension(
  HSMM,
  norm_method = 'log',
  max_components = 2,
  reduction_method = randomized_pca_monocle2, #'ICA',
  verbose = FALSE)

# 4/5 order cells to generate pseudotime
HSMM_ica <- orderCells(HSMM_ica)
HSMM_rpca <- orderCells(HSMM_rpca)

# 5/5 visuliazation
p_pseudotime_ica <- plot_cell_trajectory(
  HSMM_ica, color_by = 'Pseudotime',
  show_branch_points = FALSE) + ggtitle('ICA')
p_pseudotime_rpca <- plot_cell_trajectory(
  HSMM_rpca, color_by = 'Pseudotime',
  show_branch_points = FALSE) + ggtitle('rPCA')

p_pseudotime <- ggarrange(p_pseudotime_ica, p_pseudotime_rpca, ncol = 2)
ggsave(p_pseudotime, filename = 'monocle_custome.png',
       width = 15, height=8)
my_genes <- row.names(subset(fData(HSMM),
                             gene_short_name %in% c("CDK1", "MEF2C", "MYH3")))

p_g_ps_ica <- plot_genes_in_pseudotime(
  HSMM_ica[my_genes, ], color_by = "Hours") + ggtitle('ICA')
p_g_ps_rpca <- plot_genes_in_pseudotime(
  HSMM_rpca[my_genes, ], color_by = "Hours") + ggtitle('rPCA')
p_g_ps <- ggarrange(p_g_ps_ica, p_g_ps_rpca, ncol=2)
ggsave(p_g_ps, filename = 'monocle_custome_genes.png',
       width = 10, height = 10)