Diffusion Map Projection
Fiona Hamey
8 June 2016
This is a demonstration of how to project new data sets onto the diffusion map of adult haematopoiesis. This projection is performed in the R programming environment. You will need the destiny package installed. This can be downloaded from Bioconductor.
library(destiny)The next step is to load the data for the 1,656 blood stem and progenitor cells. These can be downloaded from http://blood.stemcells.cam.ac.uk/single_cell_atlas.html. Rows should be labelled with Ensembl gene IDs and columns with cell names.
originalExpressionMatrix <- read.table('normalisedCountsVariableGenes.txt', header = T, row.names = 1)
originalExpressionMatrix[1:5,1:5]## HSPC_025 HSPC_031 HSPC_037 LT.HSC_001 HSPC_001
## ENSMUSG00000030159 0.000000 0.0000 0.0000 0.0000 0.0000
## ENSMUSG00000053470 28.683808 116.5984 0.0000 0.0000 0.0000
## ENSMUSG00000041729 1.687283 0.0000 119.5414 288.6915 710.9228
## ENSMUSG00000048489 0.000000 0.0000 0.0000 0.0000 0.0000
## ENSMUSG00000046080 0.000000 0.0000 0.0000 0.0000 0.0000We now perform the diffusion map dimensionality reduction on these data using functions from the destiny package. The log transformation acts as a variance stabilising transformation. Here for simplicity we have used the plot function to plot the diffusion maps. To alter the viewing angle or interactively rotate the 3D plot, packages such as rgl and plot3D can be used.
logOriginalExpression <- log2(originalExpressionMatrix + 1)
dm <- DiffusionMap(t(logOriginalExpression), distance = "cosine", sigma = .16)
plot(dm, c(3,2,1), pch=20, col="grey")
The diffusion map can be coloured in by features such as cluster ID.
Next we load and prepare the new data set that we wish to project. Here we are demonstrating projection of the data from (Grover et al., 2016). The normalised expression matrix for this demonstration can be downloaded from http://blood.stemcells.cam.ac.uk/single_cell_atlas.html, or you can use your own data for projection.
newExpressionMatrix <- read.table('grover_expression.txt', header = T)To perform the projection we require both new and old data to have the same set of genes. Then we can perform the projection using the dm.predict function from the destiny package.
newExpressionMatrix <- newExpressionMatrix[rownames(originalExpressionMatrix), ]
logNewExpression <- log2(newExpressionMatrix + 1)
dmProject <- dm.predict(dm, t(logNewExpression))
plot(dm, c(3,2,1), col = "grey", new.dcs=dmProject, pch=20, col.new = "red")
devtools::session_info()## Session info --------------------------------------------------------------## setting value
## version R version 3.3.0 beta (2016-04-10 r70470)
## system x86_64, darwin13.4.0
## ui X11
## language (EN)
## collate en_GB.UTF-8
## tz Europe/London
## date 2016-06-09## Packages ------------------------------------------------------------------## package * version date source
## acepack 1.3-3.3 2014-11-24 CRAN (R 3.3.0)
## Biobase * 2.31.3 2016-01-14 Bioconductor
## BiocGenerics * 0.17.5 2016-04-11 Bioconductor
## car 2.1-2 2016-03-25 CRAN (R 3.3.0)
## chron 2.3-47 2015-06-24 CRAN (R 3.3.0)
## class 7.3-14 2015-08-30 CRAN (R 3.3.0)
## cluster 2.0.4 2016-04-18 CRAN (R 3.3.0)
## colorspace 1.2-6 2015-03-11 CRAN (R 3.3.0)
## data.table 1.9.6 2015-09-19 CRAN (R 3.3.0)
## DEoptimR 1.0-4 2015-10-23 CRAN (R 3.3.0)
## destiny * 1.1.1 2016-03-31 Bioconductor
## devtools 1.11.1 2016-04-21 CRAN (R 3.3.0)
## digest 0.6.9 2016-01-08 CRAN (R 3.3.0)
## e1071 1.6-7 2015-08-05 CRAN (R 3.3.0)
## evaluate 0.9 2016-04-29 CRAN (R 3.3.0)
## FNN 1.1 2013-07-31 CRAN (R 3.3.0)
## foreign 0.8-66 2015-08-19 CRAN (R 3.3.0)
## formatR 1.4 2016-05-09 CRAN (R 3.3.0)
## Formula 1.2-1 2015-04-07 CRAN (R 3.3.0)
## ggplot2 2.1.0 2016-03-01 CRAN (R 3.3.0)
## gridExtra 2.2.1 2016-02-29 CRAN (R 3.3.0)
## gtable 0.2.0 2016-02-26 CRAN (R 3.3.0)
## Hmisc 3.17-4 2016-05-02 CRAN (R 3.3.0)
## htmltools 0.3.5 2016-03-21 CRAN (R 3.3.0)
## igraph 1.0.1 2015-06-26 CRAN (R 3.3.0)
## knitr 1.13 2016-05-09 CRAN (R 3.3.0)
## lattice 0.20-33 2015-07-14 CRAN (R 3.3.0)
## latticeExtra 0.6-28 2016-02-09 CRAN (R 3.3.0)
## lme4 1.1-12 2016-04-16 CRAN (R 3.3.0)
## lmtest 0.9-34 2015-06-06 CRAN (R 3.3.0)
## magrittr 1.5 2014-11-22 CRAN (R 3.3.0)
## MASS 7.3-45 2015-11-10 CRAN (R 3.3.0)
## Matrix 1.2-6 2016-05-02 CRAN (R 3.3.0)
## MatrixModels 0.4-1 2015-08-22 CRAN (R 3.3.0)
## memoise 1.0.0 2016-01-29 CRAN (R 3.3.0)
## mgcv 1.8-12 2016-03-03 CRAN (R 3.3.0)
## minqa 1.2.4 2014-10-09 CRAN (R 3.3.0)
## munsell 0.4.3 2016-02-13 CRAN (R 3.3.0)
## nlme 3.1-127 2016-04-16 CRAN (R 3.3.0)
## nloptr 1.0.4 2014-08-04 CRAN (R 3.3.0)
## nnet 7.3-12 2016-02-02 CRAN (R 3.3.0)
## pbkrtest 0.4-6 2016-01-27 CRAN (R 3.3.0)
## plyr 1.8.4 2016-06-08 CRAN (R 3.3.0)
## proxy 0.4-15 2015-07-08 CRAN (R 3.3.0)
## quantreg 5.21 2016-02-13 CRAN (R 3.3.0)
## RColorBrewer 1.1-2 2014-12-07 CRAN (R 3.3.0)
## Rcpp 0.12.5 2016-05-14 CRAN (R 3.3.0)
## RcppEigen 0.3.2.8.1 2016-03-01 CRAN (R 3.3.0)
## rmarkdown 0.9.6 2016-05-01 CRAN (R 3.3.0)
## robustbase 0.92-5 2015-07-22 CRAN (R 3.3.0)
## rpart 4.1-10 2015-06-29 CRAN (R 3.3.0)
## scales 0.4.0 2016-02-26 CRAN (R 3.3.0)
## scatterplot3d 0.3-36 2015-07-30 CRAN (R 3.3.0)
## sp 1.2-3 2016-04-14 CRAN (R 3.3.0)
## SparseM 1.7 2015-08-15 CRAN (R 3.3.0)
## stringi 1.0-1 2015-10-22 CRAN (R 3.3.0)
## stringr 1.0.0 2015-04-30 CRAN (R 3.3.0)
## survival 2.39-2 2016-04-16 CRAN (R 3.3.0)
## vcd 1.4-1 2015-07-08 CRAN (R 3.3.0)
## VIM 4.4.1 2015-09-15 CRAN (R 3.3.0)
## withr 1.0.1 2016-02-04 CRAN (R 3.3.0)
## yaml 2.1.13 2014-06-12 CRAN (R 3.3.0)
## zoo 1.7-13 2016-05-03 CRAN (R 3.3.0)