--- bibliography: ref.bib --- # (PART) Case studies {-} # Cross-annotating human pancreas {#pancreas-case-study} ## Loading the data We load the @muraro2016singlecell dataset as our reference, removing unlabelled cells or cells without a clear label. ```r library(scRNAseq) sceM <- MuraroPancreasData() sceM <- sceM[,!is.na(sceM$label) & sceM$label!="unclear"] ``` We compute log-expression values for use in marker detection inside `SingleR()`. ```r library(scater) sceM <- logNormCounts(sceM) ``` We examine the distribution of labels in this reference. ```r table(sceM$label) ``` ``` ## ## acinar alpha beta delta duct endothelial ## 219 812 448 193 245 21 ## epsilon mesenchymal pp ## 3 80 101 ``` We load the @grun2016denovo dataset as our test, applying some basic quality control to remove low-quality cells in some of the batches (see [here](https://osca.bioconductor.org/grun-human-pancreas-cel-seq2.html#quality-control-8) for details). ```r sceG <- GrunPancreasData() sceG <- addPerCellQC(sceG) qc <- quickPerCellQC(colData(sceG), percent_subsets="altexps_ERCC_percent", batch=sceG$donor, subset=sceG$donor %in% c("D17", "D7", "D2")) sceG <- sceG[,!qc$discard] ``` Technically speaking, the test dataset does not need log-expression values but we compute them anyway for convenience. ```r sceG <- logNormCounts(sceG) ``` ## Applying the annotation We apply `SingleR()` with Wilcoxon rank sum test-based marker detection to annotate the Grun dataset with the Muraro labels. ```r library(SingleR) pred.grun <- SingleR(test=sceG, ref=sceM, labels=sceM$label, de.method="wilcox") ``` We examine the distribution of predicted labels: ```r table(pred.grun$labels) ``` ``` ## ## acinar alpha beta delta duct endothelial ## 277 203 181 50 306 5 ## epsilon mesenchymal pp ## 1 22 19 ``` We can also examine the number of discarded cells for each label: ```r table(Label=pred.grun$labels, Lost=is.na(pred.grun$pruned.labels)) ``` ``` ## Lost ## Label FALSE TRUE ## acinar 251 26 ## alpha 198 5 ## beta 180 1 ## delta 49 1 ## duct 301 5 ## endothelial 4 1 ## epsilon 1 0 ## mesenchymal 22 0 ## pp 17 2 ``` ## Diagnostics We visualize the assignment scores for each label in Figure \@ref(fig:unref-pancreas-score-heatmap). ```r plotScoreHeatmap(pred.grun) ```
Heatmap of the (normalized) assignment scores for each cell (column) in the Grun test dataset with respect to each label (row) in the Muraro reference dataset. The final assignment for each cell is shown in the annotation bar at the top.

(\#fig:unref-pancreas-score-heatmap)Heatmap of the (normalized) assignment scores for each cell (column) in the Grun test dataset with respect to each label (row) in the Muraro reference dataset. The final assignment for each cell is shown in the annotation bar at the top.

The delta for each cell is visualized in Figure \@ref(fig:unref-pancreas-delta-dist). ```r plotDeltaDistribution(pred.grun) ```
Distributions of the deltas for each cell in the Grun dataset assigned to each label in the Muraro dataset. Each cell is represented by a point; low-quality assignments that were pruned out are colored in orange.

(\#fig:unref-pancreas-delta-dist)Distributions of the deltas for each cell in the Grun dataset assigned to each label in the Muraro dataset. Each cell is represented by a point; low-quality assignments that were pruned out are colored in orange.

Finally, we visualize the heatmaps of the marker genes for each label in Figure \@ref(fig:unref-pancreas-marker-heat). ```r library(scater) collected <- list() all.markers <- metadata(pred.grun)$de.genes sceG$labels <- pred.grun$labels for (lab in unique(pred.grun$labels)) { collected[[lab]] <- plotHeatmap(sceG, silent=TRUE, order_columns_by="labels", main=lab, features=unique(unlist(all.markers[[lab]])))[[4]] } do.call(gridExtra::grid.arrange, collected) ```
Heatmaps of log-expression values in the Grun dataset for all marker genes upregulated in each label in the Muraro reference dataset. Assigned labels for each cell are shown at the top of each plot.

(\#fig:unref-pancreas-marker-heat)Heatmaps of log-expression values in the Grun dataset for all marker genes upregulated in each label in the Muraro reference dataset. Assigned labels for each cell are shown at the top of each plot.

## Comparison to clusters For comparison, we will perform a quick unsupervised analysis of the Grun dataset. We model the variances using the spike-in data and we perform graph-based clustering (increasing the resolution by dropping `k=5`). ```r library(scran) decG <- modelGeneVarWithSpikes(sceG, "ERCC") set.seed(1000100) sceG <- denoisePCA(sceG, decG) library(bluster) sceG$cluster <- clusterRows(reducedDim(sceG), NNGraphParam(k=5)) ``` We see that the clusters map reasonably well to the labels in Figure \@ref(fig:unref-pancreas-label-clusters). ```r tab <- table(cluster=sceG$cluster, label=pred.grun$labels) pheatmap::pheatmap(log10(tab+10)) ```
Heatmap of the log-transformed number of cells in each combination of label (column) and cluster (row) in the Grun dataset.

(\#fig:unref-pancreas-label-clusters)Heatmap of the log-transformed number of cells in each combination of label (column) and cluster (row) in the Grun dataset.

We proceed to the most important part of the analysis. Yes, that's right, the $t$-SNE plot (Figure \@ref(fig:unref-pancreas-label-tsne)). ```r set.seed(101010100) sceG <- runTSNE(sceG, dimred="PCA") plotTSNE(sceG, colour_by="cluster", text_colour="red", text_by=I(pred.grun$labels)) ```
$t$-SNE plot of the Grun dataset, where each point is a cell and is colored by the assigned cluster. Reference labels from the Muraro dataset are also placed on the median coordinate across all cells assigned with that label.

(\#fig:unref-pancreas-label-tsne)$t$-SNE plot of the Grun dataset, where each point is a cell and is colored by the assigned cluster. Reference labels from the Muraro dataset are also placed on the median coordinate across all cells assigned with that label.

## Session information {-}
``` R version 4.0.4 (2021-02-15) Platform: x86_64-pc-linux-gnu (64-bit) Running under: Ubuntu 20.04.2 LTS Matrix products: default BLAS: /home/biocbuild/bbs-3.12-books/R/lib/libRblas.so LAPACK: /home/biocbuild/bbs-3.12-books/R/lib/libRlapack.so locale: [1] LC_CTYPE=en_US.UTF-8 LC_NUMERIC=C [3] LC_TIME=en_US.UTF-8 LC_COLLATE=C [5] LC_MONETARY=en_US.UTF-8 LC_MESSAGES=en_US.UTF-8 [7] LC_PAPER=en_US.UTF-8 LC_NAME=C [9] LC_ADDRESS=C LC_TELEPHONE=C [11] LC_MEASUREMENT=en_US.UTF-8 LC_IDENTIFICATION=C attached base packages: [1] parallel stats4 stats graphics grDevices utils datasets [8] methods base other attached packages: [1] bluster_1.0.0 scran_1.18.5 [3] SingleR_1.4.1 scater_1.18.6 [5] ggplot2_3.3.3 scRNAseq_2.4.0 [7] SingleCellExperiment_1.12.0 SummarizedExperiment_1.20.0 [9] Biobase_2.50.0 GenomicRanges_1.42.0 [11] GenomeInfoDb_1.26.4 IRanges_2.24.1 [13] S4Vectors_0.28.1 BiocGenerics_0.36.0 [15] MatrixGenerics_1.2.1 matrixStats_0.58.0 [17] BiocStyle_2.18.1 rebook_1.0.0 loaded via a namespace (and not attached): [1] AnnotationHub_2.22.0 BiocFileCache_1.14.0 [3] igraph_1.2.6 lazyeval_0.2.2 [5] BiocParallel_1.24.1 digest_0.6.27 [7] ensembldb_2.14.0 htmltools_0.5.1.1 [9] viridis_0.5.1 fansi_0.4.2 [11] magrittr_2.0.1 memoise_2.0.0 [13] limma_3.46.0 Biostrings_2.58.0 [15] askpass_1.1 prettyunits_1.1.1 [17] colorspace_2.0-0 blob_1.2.1 [19] rappdirs_0.3.3 xfun_0.22 [21] dplyr_1.0.5 callr_3.5.1 [23] crayon_1.4.1 RCurl_1.98-1.3 [25] jsonlite_1.7.2 graph_1.68.0 [27] glue_1.4.2 gtable_0.3.0 [29] zlibbioc_1.36.0 XVector_0.30.0 [31] DelayedArray_0.16.2 BiocSingular_1.6.0 [33] scales_1.1.1 pheatmap_1.0.12 [35] DBI_1.1.1 edgeR_3.32.1 [37] Rcpp_1.0.6 viridisLite_0.3.0 [39] xtable_1.8-4 progress_1.2.2 [41] dqrng_0.2.1 bit_4.0.4 [43] rsvd_1.0.3 httr_1.4.2 [45] RColorBrewer_1.1-2 ellipsis_0.3.1 [47] farver_2.1.0 pkgconfig_2.0.3 [49] XML_3.99-0.6 scuttle_1.0.4 [51] CodeDepends_0.6.5 sass_0.3.1 [53] dbplyr_2.1.0 locfit_1.5-9.4 [55] utf8_1.2.1 labeling_0.4.2 [57] tidyselect_1.1.0 rlang_0.4.10 [59] later_1.1.0.1 AnnotationDbi_1.52.0 [61] munsell_0.5.0 BiocVersion_3.12.0 [63] tools_4.0.4 cachem_1.0.4 [65] generics_0.1.0 RSQLite_2.2.4 [67] ExperimentHub_1.16.0 evaluate_0.14 [69] stringr_1.4.0 fastmap_1.1.0 [71] yaml_2.2.1 processx_3.4.5 [73] knitr_1.31 bit64_4.0.5 [75] purrr_0.3.4 AnnotationFilter_1.14.0 [77] sparseMatrixStats_1.2.1 mime_0.10 [79] xml2_1.3.2 biomaRt_2.46.3 [81] compiler_4.0.4 beeswarm_0.3.1 [83] curl_4.3 interactiveDisplayBase_1.28.0 [85] statmod_1.4.35 tibble_3.1.0 [87] bslib_0.2.4 stringi_1.5.3 [89] highr_0.8 ps_1.6.0 [91] GenomicFeatures_1.42.2 lattice_0.20-41 [93] ProtGenerics_1.22.0 Matrix_1.3-2 [95] vctrs_0.3.6 pillar_1.5.1 [97] lifecycle_1.0.0 BiocManager_1.30.10 [99] jquerylib_0.1.3 BiocNeighbors_1.8.2 [101] cowplot_1.1.1 bitops_1.0-6 [103] irlba_2.3.3 httpuv_1.5.5 [105] rtracklayer_1.50.0 R6_2.5.0 [107] bookdown_0.21 promises_1.2.0.1 [109] gridExtra_2.3 vipor_0.4.5 [111] codetools_0.2-18 assertthat_0.2.1 [113] openssl_1.4.3 withr_2.4.1 [115] GenomicAlignments_1.26.0 Rsamtools_2.6.0 [117] GenomeInfoDbData_1.2.4 hms_1.0.0 [119] grid_4.0.4 beachmat_2.6.4 [121] rmarkdown_2.7 DelayedMatrixStats_1.12.3 [123] Rtsne_0.15 shiny_1.6.0 [125] ggbeeswarm_0.6.0 ```