library(heatmaply)
library(dendextend)
Here we wish to reproduce figure 3: “Figure 3: Heatmap analysis of the most abundant bacterial genera detected across all samples.”.
The heatmap depicts the relative percentage of 16S rRNA gene sequences assigned to each bacterial genus (y axis) across the 8 samples analysed (x axis). The heatmap colors represent the relative percentage of the microbial genus assignments within each sample. Square colors shifted towards bright blue indicate higher abundance. The relative abundance values of each genus for each sample are reported in Suppl. Tab. 4.
The data was shared by Alfano, Niccolo. It shows for each genus of bacteria detected the relative abundance in the different samples used (eye, mouth, rectum, faeces of 2 koalas). Please keep in mind that the sum of the abundances of the genera could be less than 100%, cause the genera with an abundance <0.1% of the total read count were removed from the table to simplify the visualization of the results.
# koala_gut_genes <- read.csv("Figure 3_table.csv")
# rownames(koala_gut_genes) <- koala_gut_genes[,1]
# koala_gut_genes <- koala_gut_genes[,-1]
# http://r-pkgs.had.co.nz/data.html#data-data
# devtools::use_data(koala_gut_genes)
library(heatmaplyExamples)
data(koala_gut_genes)
head(koala_gut_genes)
#> Eye_SN265 Eye_SN241 Mouth_SN265 Mouth_SN241
#> Corynebacterium 14.84 10.46 0.00 0.91
#> uncultured_Corynebacteriaceae 0.06 6.58 0.00 0.06
#> Propionibacterium 0.00 0.19 0.00 6.85
#> uncultured_BD2-2 0.13 0.00 0.00 0.00
#> Bacteroides 3.16 6.75 0.03 14.36
#> Barnesiella 0.31 0.16 0.00 0.32
#> Rectum_SN265 Rectum_SN241 Faeces_SN265
#> Corynebacterium 0.10 0.72 0.00
#> uncultured_Corynebacteriaceae 0.06 0.43 0.00
#> Propionibacterium 6.08 0.35 0.00
#> uncultured_BD2-2 3.38 0.00 1.91
#> Bacteroides 52.04 0.29 71.87
#> Barnesiella 4.13 0.01 3.99
#> Faeces_SN241
#> Corynebacterium 0.00
#> uncultured_Corynebacteriaceae 0.00
#> Propionibacterium 0.00
#> uncultured_BD2-2 0.00
#> Bacteroides 23.14
#> Barnesiella 0.74
Note that the color palette used here is not linear in the data.
library(heatmaply)
library(RColorBrewer)
library(scales) # for using a scaling of the colors which is similar to what appeared in the paper.
# display.brewer.all()
black_blue_fun <- gradient_n_pal(c("black", "royalblue1"), values = rescale(c(0,1:5,10,15,20,30,50,60,70,80,100)))
black_blue <- black_blue_fun(seq(0,1,length.out = 100))
heatmaply(koala_gut_genes, dendrogram=FALSE,
col = black_blue, limits = c(0,100))
We can have the same heatmap but using the viridis color palette:
library(viridis)
viridis_fun <- gradient_n_pal(viridis(100), values = rescale(c(0,1:5,10,15,20,30,50,60,70,80,100)))
viridis_cols <- viridis_fun(seq(0,1,length.out = 100))
heatmaply(koala_gut_genes, dendrogram=FALSE,
col = viridis_cols, limits = c(0,100))
And the same can be done, but with adding dendrograms and ordering them:
heatmaply(koala_gut_genes,
col = viridis_cols, limits = c(0,100))
sessionInfo()
#> R version 3.4.1 (2017-06-30)
#> Platform: x86_64-w64-mingw32/x64 (64-bit)
#> Running under: Windows 7 x64 (build 7601) Service Pack 1
#>
#> Matrix products: default
#>
#> locale:
#> [1] LC_COLLATE=Hebrew_Israel.1255 LC_CTYPE=Hebrew_Israel.1255
#> [3] LC_MONETARY=Hebrew_Israel.1255 LC_NUMERIC=C
#> [5] LC_TIME=Hebrew_Israel.1255
#>
#> attached base packages:
#> [1] parallel stats graphics grDevices datasets utils methods
#> [8] base
#>
#> other attached packages:
#> [1] scales_0.5.0 RColorBrewer_1.1-2
#> [3] limma_3.22.7 ALL_1.7.1
#> [5] Biobase_2.26.0 BiocGenerics_0.12.1
#> [7] dendextend_1.6.0 glmnet_2.0-10
#> [9] foreach_1.4.3 Matrix_1.2-10
#> [11] bindrcpp_0.2 heatmaplyExamples_0.2.0
#> [13] knitr_1.16 gplots_3.0.1
#> [15] heatmaply_0.11.0 viridis_0.4.0
#> [17] viridisLite_0.2.0 plotly_4.7.1
#> [19] ggplot2_2.2.1.9000 installr_0.19.0
#> [21] stringr_1.2.0
#>
#> loaded via a namespace (and not attached):
#> [1] httr_1.2.1 tidyr_0.6.3 jsonlite_1.5
#> [4] gtools_3.5.0 shiny_1.0.3 assertthat_0.2.0
#> [7] stats4_3.4.1 yaml_2.1.14 robustbase_0.92-7
#> [10] backports_1.1.0 lattice_0.20-35 glue_1.1.1
#> [13] digest_0.6.12 colorspace_1.3-2 httpuv_1.3.5
#> [16] htmltools_0.3.6 plyr_1.8.4 devtools_1.13.2
#> [19] pkgconfig_2.0.1 xtable_1.8-2 purrr_0.2.2.2
#> [22] mvtnorm_1.0-6 gdata_2.18.0 whisker_0.3-2
#> [25] tibble_1.3.4 mgcv_1.8-17 withr_2.0.0
#> [28] nnet_7.3-12 lazyeval_0.2.0 mime_0.5
#> [31] magrittr_1.5 mclust_5.3 memoise_1.1.0
#> [34] evaluate_0.10.1 nlme_3.1-131 MASS_7.3-47
#> [37] class_7.3-14 tools_3.4.1 registry_0.3
#> [40] data.table_1.10.4 trimcluster_0.1-2 kernlab_0.9-25
#> [43] munsell_0.4.3 cluster_2.0.6 fpc_2.1-10
#> [46] compiler_3.4.1 caTools_1.17.1 rlang_0.1.2
#> [49] grid_3.4.1 iterators_1.0.8 htmlwidgets_0.9
#> [52] crosstalk_1.0.0 labeling_0.3 bitops_1.0-6
#> [55] rmarkdown_1.6 gtable_0.2.0 codetools_0.2-15
#> [58] flexmix_2.3-14 reshape2_1.4.2 TSP_1.1-5
#> [61] R6_2.2.2 seriation_1.2-2 gridExtra_2.2.1
#> [64] prabclus_2.2-6 dplyr_0.7.3.9000 bindr_0.1
#> [67] rprojroot_1.2 KernSmooth_2.23-15 modeltools_0.2-21
#> [70] stringi_1.1.5 Rcpp_0.12.12 gclus_1.3.1
#> [73] DEoptimR_1.0-8 diptest_0.75-7