<jupyter_start><jupyter_text>*First compiled: September 17, 2017.*  Monocle 2 for data of [Paul *et al.*, Cell (2015)](https://doi.org/10.1016/j.cell.2015.11.013): alternative preprocessing Here we use an alternative preprocessing as in [*monocle2_original*](monocle2_original.ipynb), which reproduces  the notebook that accompanies the Monocle 2 publication of [Qiu *et al.*, Nat. Meth. (2017)](https://doi.org/10.1038/nmeth.4402), available from [here](https://github.com/cole-trapnell-lab/monocle2-rge-paper/blob/26ce18c97f22e9488eb85971a23311f5252bfce5/Paul_dataset_analysis_final.ipynb).Instead of manually removing a cluster of lympoid cells, we keep this cluster. Everything else is exactly as in the original notebook.<jupyter_code>rm(list = ls()) # clear the environment 
options(warn=-1) # turn off warning message globally 
suppressMessages(library(monocle))
suppressMessages(library(plyr))
suppressMessages(library(dplyr))<jupyter_output><empty_output><jupyter_text>Load the data This RData is from <NAME> in the Theis lab https://github.com/theislab/scAnalysisTutorial. It's the same data as the one available within Scanpy via [*sc.examples.paul15()*](https://scanpy.readthedocs.io/en/latest/api/scanpy.api.datasets.paul15.html).<jupyter_code>load('./data/Paul_Cell_MARSseq_GSE72857.RData') 
# the following code is used to select feature genes used by Maren 
gene.names <-sapply(strsplit(rownames(data.debatched), ";"), "[", 1)
is.informative <- gene.names %in% info.genes[order(info.genes)]
data.info.genes <- data.debatched[is.informative,]
rownames(data.info.genes) <- gene.names[is.informative]<jupyter_output><empty_output><jupyter_text>Create the cell datasets In constrast to the preprocessing as done for Supplemental Figure 16 of [Qiu *et al.*, Nat. Meth. (2017)](https://doi.org/10.1038/nmeth.4402) and as reproduced here in [*monocle2_original*](monocle2_original.ipynb), now, we do not filter out the cluster of lymphoid cells.<jupyter_code>previous_time <- proc.time()[3]
################################################################################################################################################
# obtain this mat file from Ido Amit group  
MAP_cells_clusters <- read.csv('./data/MAP.csv', header = F)
row.names(MAP_cells_clusters) <- MAP_cells_clusters$V1

#filtering cells to include only the ones which were assigned a cluster id: 
valid_subset_GSE72857_exprs <- read.table('./data/GSE72857_umitab.txt', header = T, row.names = 1)
design_mat <- read.table('./data/GSE72857_experimental_design.txt', header = T, row.names = 1, skip = 19, sep = '\t')
design_mat$cluster <- MAP_cells_clusters[row.names(design_mat), 'V2']
valid_design_mat <- subset(design_mat, !is.na(cluster))

# Get the intersect gene used by <NAME> and the genes we have 
common_genes <- rownames(valid_subset_GSE72857_exprs)[rownames(valid_subset_GSE72857_exprs) %in% info.genes]
fd <- new("AnnotatedDataFrame", data = data.frame(gene_short_name = common_genes, row.names = common_genes))
pd <- new("AnnotatedDataFrame", data = valid_design_mat)

# create a CDS with data.info.genes 
valid_subset_GSE72857_cds <- newCellDataSet(as(as.matrix(data.info.genes[common_genes, ]), 'sparseMatrix'), 
                                            phenoData = pd, 
                                            featureData = fd,
                                            lowerDetectionLimit=1,
                                            expressionFamily=negbinomial.size())
valid_subset_GSE72857_cds <- estimateSizeFactors(valid_subset_GSE72857_cds)
valid_subset_GSE72857_cds <- estimateDispersions(valid_subset_GSE72857_cds)

pData(valid_subset_GSE72857_cds)$cell_type <- revalue(as.character(pData(valid_subset_GSE72857_cds)$cluster), 
                                                      c("1" = 'erythroid', "2" = 'erythroid', "3" = 'erythroid', "4" = 'erythroid', "5" = 'erythroid', "6" = 'erythroid', 
                                                        "7" = 'CMP', "8" = 'CMP', "9" = 'CMP', "10" = 'CMP',
                                                        "11" = 'DC', 
                                                        "12" = 'GMP', "13" = 'GMP', "14" = 'GMP', "15" = 'GMP', "16" = 'GMP', "17" = 'GMP', "18" = 'GMP', 
                                                        "19" = 'lymphoid'))

# remove all lymphoid cells as they do not belong to myeloid lineage
# HERE, we commented out the next line; we do NOT remove the lymphoid cells
# valid_subset_GSE72857_cds <- valid_subset_GSE72857_cds[, pData(valid_subset_GSE72857_cds)$cell_type != 'lymphoid']
proc.time()[3] - previous_time<jupyter_output>Removing 23 outliers<jupyter_text>Reconstructing the trajectory with Monocle 2  The variance explained by each component in PCA.<jupyter_code>options(repr.plot.width=4, repr.plot.height=3)
plot_pc_variance_explained(valid_subset_GSE72857_cds) + geom_vline(xintercept = 10)<jupyter_output>Warning message in (function (A, nv = 5, nu = nv, maxit = 100, work = nv + 7, reorth = TRUE, :
“did not converge--results might be invlaid!; try increasing maxit or fastpath=FALSE”<jupyter_text>The actual computation.<jupyter_code>previous_time <- proc.time()[3]
valid_subset_GSE72857_cds2 <- reduceDimension(valid_subset_GSE72857_cds,
                                              norm_method = 'log', verbose = F, max_components = 10) 
valid_subset_GSE72857_cds2 <- orderCells(valid_subset_GSE72857_cds2, reverse = T)
proc.time()[3] - previous_time<jupyter_output><empty_output><jupyter_text>Similar analysis as for Supplemental Figure 16 of [Qiu *et al.*, Nat. Meth. (2017)](https://doi.org/10.1038/nmeth.4402)  Analogous to Supplemental Figure 16a Without removing the lymphoid cells, the figure fundamentally differs from the original [notebook](https://github.com/cole-trapnell-lab/monocle2-rge-paper/blob/26ce18c97f22e9488eb85971a23311f5252bfce5/Paul_dataset_analysis_final.ipynb).<jupyter_code>detailed_cell_type_color <- c("B" = "#E088B8", "DC" = "#46C7EF", "Eos" = "#EFAD1E", "Ery" = "#8CB3DF", "Mo" = "#53C0AD", "MEP" = "#4EB859", "GMP" = "#D097C4", "MK" = "#ACC436", "Neu" = "#F5918A", "lymphoid" = "#FF0000")

pData(valid_subset_GSE72857_cds2)$cell_type2 <- revalue(as.character(pData(valid_subset_GSE72857_cds2)$cluster), 
                                                        c("1" = 'Ery', "2" = 'Ery', "3" = 'Ery', "4" = 'Ery', "5" = 'Ery', "6" = 'Ery', 
                                                          "7" = 'MEP', "8" = 'MK', "9" = 'GMP', "10" = 'GMP',
                                                          "11" = 'DC', 
                                                          "12" = 'B', "13" = 'B', "14" = 'Mo', "15" = 'Mo', "16" = 'Neu', "17" = 'Neu', "18" = 'Eos', 
                                                          "19" = 'lymphoid'))

options(repr.plot.width=3, repr.plot.height=4)
plot_complex_cell_trajectory(valid_subset_GSE72857_cds2, color_by = 'State', show_branch_points = T, 
                             cell_size = 0.5, cell_link_size = 0.3, root_states = c(15))<jupyter_output><empty_output><jupyter_text>Analogous to Supplemental Figure 16b Without removing the lymphoid cells, the figure fundamentally differs from the original [notebook](https://github.com/cole-trapnell-lab/monocle2-rge-paper/blob/26ce18c97f22e9488eb85971a23311f5252bfce5/Paul_dataset_analysis_final.ipynb).<jupyter_code>options(repr.plot.width=4, repr.plot.height=3)
plot_complex_cell_trajectory(valid_subset_GSE72857_cds2, color_by = 'as.factor(cell_type2)', show_branch_points = T, 
                             cell_size = 0.5, cell_link_size = 0.3, root_states = c(15)) + scale_size(range = c(0.2, 0.2)) +
  theme(axis.text.x = element_text(angle = 30, hjust = 1)) +
  theme (legend.position="left", legend.title=element_blank()) + scale_color_manual(values = detailed_cell_type_color)<jupyter_output><empty_output><jupyter_text>Session Info<jupyter_code>sessionInfo()<jupyter_output><empty_output>
{'fraction_non_alphanumeric': 0.1473247927656368, 'fraction_numerical': 0.055639286611404166, 'mean_word_length': 3.4861971830985916, 'pattern_counts': {'":': 0, '<': 57, '<?xml version=': 0, '>': 33, 'https://': 9, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 30, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '8', 'dir': 'jupyter-structured-clean-dedup', 'id': '221037', 'license': 'BSD-3-Clause', 'n_tokens_mistral': 2982, 'n_tokens_neox': 2585, 'n_words': 541, 'path': 'comparisons/paul15_monocle2/monocle2_alternative.ipynb', 'repo_name': 'evanbiederstedt/paga', 'size': '7976'}