- download: null
- edit: https://github.com/OSCA-source/OSCA.basic/edit/devel/inst/book/%s
- view: https://github.com/OSCA-source/OSCA.basic/blob/devel/inst/book/%s
- history: https://github.com/OSCA-source/OSCA.basic/commits/devel/inst/book/%s
- sharing:
- github: yes
- facebook: yes
- twitter: yes
- linkedin: yes
diff --git a/inst/book/cell-annotation.Rmd b/inst/book/cell-annotation.Rmd
deleted file mode 100644
index 5090cf4..0000000
--- a/inst/book/cell-annotation.Rmd
+++ /dev/null
@@ -1,454 +0,0 @@
----
-output: html_document
-bibliography: ref.bib
----
-
-# Cell type annotation
-
-```{r setup, echo=FALSE, results="asis"}
-library(rebook)
-chapterPreamble()
-```
-
-## Motivation
-
-The most challenging task in scRNA-seq data analysis is arguably the interpretation of the results.
-Obtaining clusters of cells is fairly straightforward, but it is more difficult to determine what biological state is represented by each of those clusters.
-Doing so requires us to bridge the gap between the current dataset and prior biological knowledge, and the latter is not always available in a consistent and quantitative manner.
-Indeed, even the concept of a "cell type" is [not clearly defined](https://doi.org/10.1016/j.cels.2017.03.006), with most practitioners possessing a "I'll know it when I see it" intuition that is not amenable to computational analysis.
-As such, interpretation of scRNA-seq data is often manual and a common bottleneck in the analysis workflow.
-
-To expedite this step, we can use various computational approaches that exploit prior information to assign meaning to an uncharacterized scRNA-seq dataset.
-The most obvious sources of prior information are the curated gene sets associated with particular biological processes, e.g., from the Gene Ontology (GO) or the Kyoto Encyclopedia of Genes and Genomes (KEGG) collections.
-Alternatively, we can directly compare our expression profiles to published reference datasets where each sample or cell has already been annotated with its putative biological state by domain experts.
-Here, we will demonstrate both approaches with several different scRNA-seq datasets.
-
-## Assigning cell labels from reference data
-
-### Overview
-
-A conceptually straightforward annotation approach is to compare the single-cell expression profiles with previously annotated reference datasets.
-Labels can then be assigned to each cell in our uncharacterized test dataset based on the most similar reference sample(s), for some definition of "similar".
-This is a standard classification challenge that can be tackled by standard machine learning techniques such as random forests and support vector machines.
-Any published and labelled RNA-seq dataset (bulk or single-cell) can be used as a reference, though its reliability depends greatly on the expertise of the original authors who assigned the labels in the first place.
-
-In this section, we will demonstrate the use of the `r Biocpkg("SingleR")` method [@aran2019reference] for cell type annotation.
-This method assigns labels to cells based on the reference samples with the highest Spearman rank correlations, using only the marker genes between pairs of labels to focus on the relevant differences between cell types.
-It also performs a fine-tuning step for each cell where the correlations are recomputed with just the marker genes for the top-scoring labels.
-This aims to resolve any ambiguity between those labels by removing noise from irrelevant markers for other labels.
-Further details can be found in the [_SingleR_ book](https://bioconductor.org/books/release/SingleRBook/) from which most of the examples here are derived.
-
-### Using existing references
-
-For demonstration purposes, we will use one of the 10X PBMC datasets as our test.
-While we have already applied quality control, normalization and clustering for this dataset, this is not strictly necessary.
-It is entirely possible to run `SingleR()` on the raw counts without any _a priori_ quality control
-and filter on the annotation results at one's leisure - see the book for an explanation.
-
-```{r, echo=FALSE, results="asis"}
-extractFromPackage("tenx-unfiltered-pbmc4k.Rmd", package="OSCA.workflows",
- chunk="clustering", objects="sce.pbmc")
-```
-
-```{r}
-sce.pbmc
-```
-
-The `r Biocpkg("celldex")` contains a number of curated reference datasets, mostly assembled from bulk RNA-seq or microarray data of sorted cell types.
-These references are often good enough for most applications provided that they contain the cell types that are expected in the test population.
-Here, we will use a reference constructed from Blueprint and ENCODE data [@martens2013blueprint;@encode2012integrated];
-this is obtained by calling the `BlueprintEncode()` function to construct a `SummarizedExperiment` containing log-expression values with curated labels for each sample.
-
-```{r}
-library(celldex)
-ref <- BlueprintEncodeData()
-ref
-```
-
-We call the `SingleR()` function to annotate each of our PBMCs with the main cell type labels from the Blueprint/ENCODE reference.
-This returns a `DataFrame` where each row corresponds to a cell in the test dataset and contains its label assignments.
-Alternatively, we could use the labels in `ref$label.fine`, which provide more resolution at the cost of speed and increased ambiguity in the assignments.
-
-```{r}
-library(SingleR)
-pred <- SingleR(test=sce.pbmc, ref=ref, labels=ref$label.main)
-table(pred$labels)
-```
-
-We inspect the results using a heatmap of the per-cell and label scores (Figure \@ref(fig:singler-heat-pbmc)).
-Ideally, each cell should exhibit a high score in one label relative to all of the others, indicating that the assignment to that label was unambiguous.
-This is largely the case for monocytes and B cells, whereas we see more ambiguity between CD4^+^ and CD8^+^ T cells (and to a lesser extent, NK cells).
-
-```{r singler-heat-pbmc, fig.wide=TRUE, fig.cap="Heatmap of the assignment score for each cell (column) and label (row). Scores are shown before any fine-tuning and are normalized to [0, 1] within each cell."}
-plotScoreHeatmap(pred)
-```
-
-We compare the assignments with the clustering results to determine the identity of each cluster.
-Here, several clusters are nested within the monocyte and B cell labels (Figure \@ref(fig:singler-cluster)), indicating that the clustering represents finer subdivisions within the cell types.
-Interestingly, our clustering does not effectively distinguish between CD4^+^ and CD8^+^ T cell labels.
-This is probably due to the presence of other factors of heterogeneity within the T cell subpopulation (e.g., activation) that have a stronger influence on unsupervised methods than the _a priori_ expected CD4^+^/CD8^+^ distinction.
-
-```{r singler-cluster, fig.cap="Heatmap of the distribution of cells across labels and clusters in the 10X PBMC dataset. Color scale is reported in the log~10~-number of cells for each cluster-label combination."}
-library(SingleCellExperiment)
-tab <- table(Assigned=pred$pruned.labels, Cluster=colLabels(sce.pbmc))
-
-# Adding a pseudo-count of 10 to avoid strong color jumps with just 1 cell.
-library(pheatmap)
-pheatmap(log2(tab+10), color=colorRampPalette(c("white", "blue"))(101))
-```
-
-```{r, echo=FALSE}
-xtab <- t(t(tab)/colSums(tab))*100
-# Checking that no cluster is enriched for CD4 T cells, supporting the text.
-stopifnot(all(xtab["CD4+ T-cells",] < 80, na.rm=TRUE))
-stopifnot(any(xtab["CD8+ T-cells",] > 80, na.rm=TRUE))
-```
-
-This episode highlights some of the differences between reference-based annotation and unsupervised clustering.
-The former explicitly focuses on aspects of the data that are known to be interesting, simplifying the process of biological interpretation.
-However, the cost is that the downstream analysis is restricted by the diversity and resolution of the available labels, a problem that is largely avoided by _de novo_ identification of clusters.
-We suggest applying both strategies to examine the agreement (or lack thereof) between reference label and cluster assignments.
-Any inconsistencies are not necessarily problematic due to the conceptual differences between the two approaches;
-indeed, one could use those discrepancies as the basis for further investigation to discover novel factors of variation in the data.
-
-### Using custom references
-
-We can also apply `r Biocpkg("SingleR")` to single-cell reference datasets that are curated and supplied by the user.
-This is most obviously useful when we have an existing dataset that was previously (manually) annotated
-and we want to use that knowledge to annotate a new dataset in an automated manner.
-To illustrate, we will use the @muraro2016singlecell human pancreas dataset as our reference.
-
-```{r, echo=FALSE, results="asis"}
-extractFromPackage("muraro-pancreas.Rmd", package="OSCA.workflows",
- chunk="normalization", objects="sce.muraro")
-```
-
-```{r}
-sce.muraro
-
-# Pruning out unknown or unclear labels.
-sce.muraro <- sce.muraro[,!is.na(sce.muraro$label) &
- sce.muraro$label!="unclear"]
-table(sce.muraro$label)
-```
-
-Our aim is to assign labels to our test dataset from @segerstolpe2016singlecell.
-We use the same call to `SingleR()` but with `de.method="wilcox"` to identify markers via pairwise Wilcoxon ranked sum tests between labels in the reference Muraro dataset.
-This re-uses the same machinery from Chapter \@ref(marker-detection); further options to fine-tune the test procedure can be passed via the `de.args` argument.
-
-```{r, echo=FALSE, results="asis"}
-extractFromPackage("segerstolpe-pancreas.Rmd", package="OSCA.workflows",
- chunk="normalization", objects="sce.seger")
-```
-
-```{r}
-# Converting to FPKM for a more like-for-like comparison to UMI counts.
-# However, results are often still good even when this step is skipped.
-library(AnnotationHub)
-hs.db <- AnnotationHub()[["AH73881"]]
-hs.exons <- exonsBy(hs.db, by="gene")
-hs.exons <- reduce(hs.exons)
-hs.len <- sum(width(hs.exons))
-
-library(scuttle)
-available <- intersect(rownames(sce.seger), names(hs.len))
-fpkm.seger <- calculateFPKM(sce.seger[available,], hs.len[available])
-
-pred.seger <- SingleR(test=fpkm.seger, ref=sce.muraro,
- labels=sce.muraro$label, de.method="wilcox")
-table(pred.seger$labels)
-```
-
-As it so happens, we are in the fortunate position where our test dataset also contains independently defined labels.
-We see strong consistency between the two sets of labels (Figure \@ref(fig:singler-comp-pancreas)), indicating that our automatic annotation is comparable to that generated manually by domain experts.
-
-```{r singler-comp-pancreas, fig.cap="Heatmap of the confusion matrix between the predicted labels (rows) and the independently defined labels (columns) in the Segerstolpe dataset. The color is proportinal to the log-transformed number of cells with a given combination of labels from each set."}
-tab <- table(pred.seger$pruned.labels, sce.seger$CellType)
-library(pheatmap)
-pheatmap(log2(tab+10), color=colorRampPalette(c("white", "blue"))(101))
-```
-
-```{r, echo=FALSE}
-# Checking that I'm not just talking shit.
-library(bluster)
-rand <- pairwiseRand(pred.seger$labels, sce.seger$CellType, mode="index")
-stopifnot(rand > 0.9)
-```
-
-An interesting question is - given a single-cell reference dataset, is it better to use it directly or convert it to pseudo-bulk values?
-A single-cell reference preserves the "shape" of the subpopulation in high-dimensional expression space, potentially yielding more accurate predictions when the differences between labels are subtle (or at least capturing ambiguity more accurately to avoid grossly incorrect predictions).
-However, it also requires more computational work to assign each cell in the test dataset.
-We refer to the [other book](https://ltla.github.io/SingleRBook/using-single-cell-references.html#pseudo-bulk-aggregation) for more details on how to achieve a compromise between these two concerns.
-
-## Assigning cell labels from gene sets
-
-A related strategy is to explicitly identify sets of marker genes that are highly expressed in each individual cell.
-This does not require matching of individual cells to the expression values of the reference dataset, which is faster and more convenient when only the identities of the markers are available.
-We demonstrate this approach using neuronal cell type markers derived from the @zeisel2015brain study.
-
-```{r, echo=FALSE, results="asis"}
-extractFromPackage("zeisel-brain.Rmd", package="OSCA.workflows",
- chunk="normalization", objects="sce.zeisel")
-```
-
-```{r}
-library(scran)
-wilcox.z <- pairwiseWilcox(sce.zeisel, sce.zeisel$level1class,
- lfc=1, direction="up")
-markers.z <- getTopMarkers(wilcox.z$statistics, wilcox.z$pairs,
- pairwise=FALSE, n=50)
-lengths(markers.z)
-```
-
-Our test dataset will be another brain scRNA-seq experiment from @tasic2016adult.
-
-```{r}
-library(scRNAseq)
-sce.tasic <- TasicBrainData()
-sce.tasic
-```
-
-We use the `r Biocpkg("AUCell")` package to identify marker sets that are highly expressed in each cell.
-This method ranks genes by their expression values within each cell and constructs a response curve of the number of genes from each marker set that are present with increasing rank.
-It then computes the area under the curve (AUC) for each marker set, quantifying the enrichment of those markers among the most highly expressed genes in that cell.
-This is roughly similar to performing a Wilcoxon rank sum test between genes in and outside of the set, but involving only the top ranking genes by expression in each cell.
-
-```{r}
-library(GSEABase)
-all.sets <- lapply(names(markers.z), function(x) {
- GeneSet(markers.z[[x]], setName=x)
-})
-all.sets <- GeneSetCollection(all.sets)
-
-library(AUCell)
-rankings <- AUCell_buildRankings(counts(sce.tasic),
- plotStats=FALSE, verbose=FALSE)
-cell.aucs <- AUCell_calcAUC(all.sets, rankings)
-results <- t(assay(cell.aucs))
-head(results)
-```
-
-We assign cell type identity to each cell in the test dataset by taking the marker set with the top AUC as the label for that cell.
-Our new labels mostly agree with the original annotation from @tasic2016adult, which is encouraging.
-The only exception involves misassignment of oligodendrocyte precursors to astrocytes, which may be understandable given that they are derived from a common lineage.
-In the absence of prior annotation, a more general diagnostic check is to compare the assigned labels to cluster identities, under the expectation that most cells of a single cluster would have the same label (or, if multiple labels are present, they should at least represent closely related cell states).
-
-```{r}
-new.labels <- colnames(results)[max.col(results)]
-tab <- table(new.labels, sce.tasic$broad_type)
-tab
-```
-
-```{r, echo=FALSE}
-stopifnot(tab["astrocytes_ependymal","Oligodendrocyte Precursor Cell"] > 0)
-
-# Checking that I'm not just talking shit.
-rand <- pairwiseRand(new.labels, sce.tasic$broad_type, mode="index")
-stopifnot(rand > 0.9)
-```
-
-As a diagnostic measure, we examine the distribution of AUCs across cells for each label (Figure \@ref(fig:auc-dist)).
-In heterogeneous populations, the distribution for each label should be bimodal with one high-scoring peak containing cells of that cell type and a low-scoring peak containing cells of other types.
-The gap between these two peaks can be used to derive a threshold for whether a label is "active" for a particular cell.
-(In this case, we simply take the single highest-scoring label per cell as the labels should be mutually exclusive.)
-In populations where a particular cell type is expected, lack of clear bimodality for the corresponding label may indicate that its gene set is not sufficiently informative.
-
-```{r auc-dist, results="hide", fig.width=8, fig.height=8, fig.cap="Distribution of AUCs in the Tasic brain dataset for each label in the Zeisel dataset. The blue curve represents the density estimate, the red curve represents a fitted two-component mixture of normals, the pink curve represents a fitted three-component mixture, and the grey curve represents a fitted normal distribution. Vertical lines represent threshold estimates corresponding to each estimate of the distribution."}
-par(mfrow=c(3,3))
-AUCell_exploreThresholds(cell.aucs, plotHist=TRUE, assign=TRUE)
-```
-
-Interpretation of the `r Biocpkg("AUCell")` results is most straightforward when the marker sets are mutually exclusive, as shown above for the cell type markers.
-In other applications, one might consider computing AUCs for gene sets associated with signalling or metabolic pathways.
-It is likely that multiple pathways will be active in any given cell, and it is tempting to use the AUCs to quantify this activity for comparison across cells.
-However, such comparisons must be interpreted with much caution as the AUCs are competitive values - any increase in one pathway's activity will naturally reduce the AUCs for all other pathways, potentially resulting in spurious differences across the population.
-
-As we mentioned previously, the advantage of the `r Biocpkg("AUCell")` approach is that it does not require reference expression values.
-This is particularly useful when dealing with gene sets derived from the literature or other qualitative forms of biological knowledge.
-For example, we might instead use single-cell signatures defined from MSigDB, obtained as shown below.
-
-```{r}
-# Downloading the signatures and caching them locally.
-library(BiocFileCache)
-bfc <- BiocFileCache(ask=FALSE)
-scsig.path <- bfcrpath(bfc, file.path("http://software.broadinstitute.org",
- "gsea/msigdb/supplemental/scsig.all.v1.0.symbols.gmt"))
-scsigs <- getGmt(scsig.path)
-```
-
-The flipside is that information on relative expression is lost when only the marker identities are used.
-The net effect of ignoring expression values is difficult to predict; for example, it may reduce performance for resolving more subtle cell types, but may also improve performance if the per-cell expression was too noisy to be useful.
-Performance is also highly dependent on the gene sets themselves, which may not be defined in the same context in which they are used.
-For example, applying all of the MSigDB signatures on the Muraro dataset is rather disappointing (Figure \@ref(fig:aucell-muraro-heat)), while restricting to the subset of pancreas signatures is more promising.
-
-```{r aucell-muraro-heat, fig.asp=1, fig.wide=TRUE, fig.cap="Heatmaps of the log-number of cells with each combination of known labels (columns) and assigned MSigDB signatures (rows) in the Muraro data set. The signature assigned to each cell was defined as that with the highest AUC across all (top) or all pancreas-related signatures (bottom)."}
-muraro.mat <- counts(sce.muraro)
-rownames(muraro.mat) <- rowData(sce.muraro)$symbol
-# Explicitly coerce count matrix to a dense matrix to avoid issues with
-# support for sparse matrices in AUCell.
-muraro.rankings <- AUCell_buildRankings(as.matrix(muraro.mat),
- plotStats=FALSE, verbose=FALSE)
-
-# Applying MsigDB to the Muraro dataset, because it's human:
-scsig.aucs <- AUCell_calcAUC(scsigs, muraro.rankings)
-scsig.results <- t(assay(scsig.aucs))
-full.labels <- colnames(scsig.results)[max.col(scsig.results)]
-tab <- table(full.labels, sce.muraro$label)
-fullheat <- pheatmap(log10(tab+10), color=viridis::viridis(100), silent=TRUE)
-
-# Restricting to the subset of Muraro-derived gene sets:
-scsigs.sub <- scsigs[grep("Pancreas", names(scsigs))]
-sub.aucs <- AUCell_calcAUC(scsigs.sub, muraro.rankings)
-sub.results <- t(assay(sub.aucs))
-sub.labels <- colnames(sub.results)[max.col(sub.results)]
-tab <- table(sub.labels, sce.muraro$label)
-subheat <- pheatmap(log10(tab+10), color=viridis::viridis(100), silent=TRUE)
-
-gridExtra::grid.arrange(fullheat[[4]], subheat[[4]])
-```
-
-```{r, echo=FALSE}
-# Checking that my trash talk above is correct.
-stopifnot(pairwiseRand(full.labels, sce.muraro$label, mode="index") < 0.2)
-stopifnot(pairwiseRand(sub.labels, sce.muraro$label, mode="index") > 0.8)
-```
-
-## Assigning cluster labels from markers
-
-```{r, echo=FALSE}
-chosen.text <- 2
-```
-
-Yet another strategy for annotation is to perform a gene set enrichment analysis on the marker genes defining each cluster.
-This identifies the pathways and processes that are (relatively) active in each cluster based on upregulation of the associated genes compared to other clusters.
-We demonstrate on the mouse mammary dataset from @bach2017differentiation, obtaining annotations for the marker genes that define cluster `r chosen.text`.
-Specifically, we define our marker subset as the top 100 genes with the largest median Cohen's $d$ (Chapter \@ref(marker-detection)).
-
-```{r, echo=FALSE, results="asis"}
-extractFromPackage("bach-mammary.Rmd", package="OSCA.workflows",
- chunk="clustering", objects="sce.mam")
-```
-
-```{r}
-markers.mam <- scoreMarkers(sce.mam, lfc=1)
-
-chosen <- "2"
-cur.markers <- markers.mam[[chosen]]
-is.de <- order(cur.markers$median.logFC.cohen, decreasing=TRUE)[1:100]
-cur.markers[is.de,1:4]
-```
-
-We test for enrichment of gene sets defined by the Gene Ontology (GO) project, which describe a comprehensive range of biological processes and functions.
-The simplest implementation of this approach involves calling the `goana()` function from the `r Biocpkg("limma")` package.
-This performs a hypergeometric test to identify GO terms that are overrepresented in our marker subset.
-
-```{r, echo=FALSE}
-options(width=100)
-```
-
-```{r}
-# goana() requires Entrez IDs, some of which map to multiple
-# symbols - hence the unique() in the call below.
-library(org.Mm.eg.db)
-entrez.ids <- mapIds(org.Mm.eg.db, keys=rownames(cur.markers),
- column="ENTREZID", keytype="SYMBOL")
-
-library(limma)
-go.out <- goana(unique(entrez.ids[is.de]), species="Mm",
- universe=unique(entrez.ids))
-
-# Only keeping biological process terms that are not overly general.
-go.out <- go.out[order(go.out$P.DE),]
-go.useful <- go.out[go.out$Ont=="BP" & go.out$N <= 200,]
-head(go.useful[,c(1,3,4)], 30)
-```
-
-We see an enrichment for genes involved in lipid storage and lipid synthesis.
-Given that this is a mammary gland experiment, we might guess that cluster `r chosen.text` contains luminal epithelial cells responsible for milk production and secretion.
-Indeed, a closer examination of the marker list indicates that this cluster upregulates milk proteins _Csn2_ and _Csn3_ (Figure \@ref(fig:violin-milk)).
-
-```{r, echo=FALSE}
-# Checking that the above statements are correct.
-stopifnot(c("GO:0019915", "GO:0019432") %in% head(rownames(go.useful), 30))
-markers.nolfc <- scoreMarkers(sce.mam, lfc=0)
-milk <- markers.nolfc[[chosen]][c("Csn2", "Csn3"),"min.logFC.cohen"]
-stopifnot(all(milk > 0))
-```
-
-```{r violin-milk, fig.asp=0.5, fig.wide=TRUE, fig.cap="Distribution of log-expression values for _Csn2_ and _Csn3_ in each cluster."}
-library(scater)
-plotExpression(sce.mam, features=c("Csn2", "Csn3"),
- x="label", colour_by="label")
-```
-
-Further inspection of interesting GO terms is achieved by extracting the relevant genes.
-This is usually desirable to confirm that the interpretation of the annotated biological process is appropriate.
-Many terms have overlapping gene sets, so a term may only be highly ranked because it shares genes with a more relevant term that represents the active pathway.
-
-```{r}
-# Extract symbols for each GO term; done once.
-tab <- select(org.Mm.eg.db, keytype="SYMBOL", keys=rownames(sce.mam), columns="GOALL")
-by.go <- split(tab[,1], tab[,2])
-
-# Identify genes associated with an interesting term.
-interesting <- unique(by.go[["GO:0006641"]])
-interesting.markers <- cur.markers[rownames(cur.markers) %in% interesting,]
-head(interesting.markers[order(-interesting.markers$median.logFC.cohen),1:4], 10)
-```
-
-Gene set testing of marker lists is a reliable approach for determining if pathways are up- or down-regulated between clusters.
-As the top marker genes are simply DEGs, we can directly apply well-established procedures for testing gene enrichment in DEG lists (see [here](https://bioconductor.org/packages/release/BiocViews.html#___GeneSetEnrichment) for relevant packages).
-This contrasts with the `r Biocpkg("AUCell")` approach where scores are not easily comparable across cells.
-The downside is that all conclusions are made relative to the other clusters, making it more difficult to determine cell identity if an "outgroup" is not present in the same study.
-
-## Computing gene set activities
-
-For the sake of completeness, we should mention that we can also quantify gene set activity on a per-cell level and test for differences in activity.
-This inverts the standard gene set testing procedure by combining information across genes first and then testing for differences afterwards.
-To avoid the pitfalls mentioned previously for the AUCs, we simply compute the average of the log-expression values across all genes in the set for each cell.
-This is less sensitive to the behavior of other genes in that cell (aside from composition biases, as discussed in Chapter \@ref(normalization)).
-
-```{r}
-aggregated <- sumCountsAcrossFeatures(sce.mam, by.go,
- exprs_values="logcounts", average=TRUE)
-dim(aggregated) # rows are gene sets, columns are cells
-aggregated[1:10,1:5]
-```
-
-We can then identify "differential gene set activity" between clusters by looking for significant differences in the per-set averages of the relevant cells.
-For example, we observe that cluster `r chosen.text` has the highest average expression for the triglyceride metabolic process GO term (Figure \@ref(fig:lipid-synth-violin)), consistent with the proposed identity of those cells.
-
-```{r lipid-synth-violin, fig.cap="Distribution of average log-normalized expression for genes involved in the triglyceride metabolic process, for all cells in each cluster of the mammary gland dataset."}
-plotColData(sce.mam, y=I(aggregated["GO:0006641",]), x="label")
-```
-
-```{r, echo=FALSE}
-maxed <- vapply(split(aggregated["GO:0006641",], colLabels(sce.mam)), mean, 0)
-stopifnot(identical(names(maxed)[which.max(maxed)], '2'))
-```
-
-The obvious disadvantage of this approach is that not all genes in the set may exhibit the same pattern of differences.
-Non-DE genes will add noise to the per-set average, "diluting" the strength of any differences compared to an analysis that focuses directly on the DE genes (Figure \@ref(fig:dbi-violin)).
-At worst, a gene set may contain subsets of DE genes that change in opposite directions, cancelling out any differences in the per-set average.
-This is not uncommon for gene sets that contain both positive and negative regulators of a particular biological process or pathway.
-
-```{r dbi-violin, fig.cap="Distribution of log-normalized expression values for _Dbi_ across all cells in each cluster of the mammary gland dataset."}
-# Choose the top-ranking gene in GO:0006641.
-plotExpression(sce.mam, "Dbi", x="label")
-```
-
-```{r, echo=FALSE}
-maxed <- vapply(split(logcounts(sce.mam)["Dbi",], colLabels(sce.mam)), mean, 0)
-stopifnot(identical(names(maxed)[which.max(maxed)], '2'))
-```
-
-We could attempt to use the per-set averages to identify gene sets of interest via differential testing across all possible sets, e.g., with `findMarkers()`.
-However, the highest ranking gene sets in this approach tend to be very small and uninteresting because - by definition - the pitfalls mentioned above are avoided when there is only one gene in the set.
-This is compounded by the fact that the log-fold changes in the per-set averages are difficult to interpret.
-For these reasons, we generally reserve the use of this gene set summary statistic for visualization rather than any real statistical analysis.
-
-## Session Info {-}
-
-```{r sessionInfo, echo=FALSE, results='asis'}
-prettySessionInfo()
-```
diff --git a/inst/book/cell-type-annotation.html b/inst/book/cell-type-annotation.html
new file mode 100644
index 0000000..66636f4
--- /dev/null
+++ b/inst/book/cell-type-annotation.html
@@ -0,0 +1,9 @@
+
+
+
+
+
+
+
+
+
diff --git a/inst/book/clustering.Rmd b/inst/book/clustering.Rmd
deleted file mode 100644
index dff44b3..0000000
--- a/inst/book/clustering.Rmd
+++ /dev/null
@@ -1,621 +0,0 @@
----
-bibliography: ref.bib
----
-
-# Clustering
-
-```{r setup, echo=FALSE, results="asis"}
-library(rebook)
-chapterPreamble()
-```
-
-## Overview
-
-Clustering is an unsupervised learning procedure that is used to empirically define groups of cells with similar expression profiles.
-Its primary purpose is to summarize complex scRNA-seq data into a digestible format for human interpretation.
-This allows us to describe population heterogeneity in terms of discrete labels that are easily understood, rather than attempting to comprehend the high-dimensional manifold on which the cells truly reside.
-After annotation based on marker genes, the clusters can be treated as proxies for more abstract biological concepts such as cell types or states.
-
-At this point, it is helpful to realize that clustering, like a microscope, is simply a tool to explore the data.
-We can zoom in and out by changing the resolution of the clustering parameters, and we can experiment with different clustering algorithms to obtain alternative perspectives of the data.
-This iterative approach is entirely permissible given that data exploration constitutes the majority of the scRNA-seq data analysis workflow.
-As such, questions about the "correctness" of the clusters or the "true" number of clusters are usually meaningless.
-We can define as many clusters as we like, with whatever algorithm we like - each clustering will represent its own partitioning of the high-dimensional expression space, and is as "real" as any other clustering.
-
-A more relevant question is "how well do the clusters approximate the cell types or states of interest?"
-Unfortunately, this is difficult to answer given the context-dependent interpretation of the underlying biology.
-Some analysts will be satisfied with resolution of the major cell types; other analysts may want resolution of subtypes; and others still may require resolution of different states (e.g., metabolic activity, stress) within those subtypes.
-Moreover, two clusterings can be highly inconsistent yet both valid, simply partitioning the cells based on different aspects of biology.
-Indeed, asking for an unqualified "best" clustering is akin to asking for the best magnification on a microscope without any context.
-
-Regardless of the exact method used, clustering is a critical step for extracting biological insights from scRNA-seq data.
-Here, we demonstrate the application of several commonly used methods with the 10X PBMC dataset.
-
-```{r, results='asis', echo=FALSE}
-extractFromPackage("tenx-unfiltered-pbmc4k.Rmd", package="OSCA.workflows",
- chunk="dimensionality-reduction", objects="sce.pbmc")
-```
-
-```{r}
-sce.pbmc
-```
-
-## Graph-based clustering {#clustering-graph}
-
-### Background
-
-Popularized by its use in `r CRANpkg("Seurat")`, graph-based clustering is a flexible and scalable technique for clustering large scRNA-seq datasets.
-We first build a graph where each node is a cell that is connected to its nearest neighbors in the high-dimensional space.
-Edges are weighted based on the similarity between the cells involved, with higher weight given to cells that are more closely related.
-We then apply algorithms to identify "communities" of cells that are more connected to cells in the same community than they are to cells of different communities.
-Each community represents a cluster that we can use for downstream interpretation.
-
-The major advantage of graph-based clustering lies in its scalability.
-It only requires a $k$-nearest neighbor search that can be done in log-linear time on average, in contrast to hierachical clustering methods with runtimes that are quadratic with respect to the number of cells.
-Graph construction avoids making strong assumptions about the shape of the clusters or the distribution of cells within each cluster, compared to other methods like $k$-means (that favor spherical clusters) or Gaussian mixture models (that require normality).
-From a practical perspective, each cell is forcibly connected to a minimum number of neighboring cells, which reduces the risk of generating many uninformative clusters consisting of one or two outlier cells.
-
-The main drawback of graph-based methods is that, after graph construction, no information is retained about relationships beyond the neighboring cells^[Sten Linarrsson talked about this in SCG2018, but I don't know where that work ended up. So this is what passes as a reference for the time being.].
-This has some practical consequences in datasets that exhibit differences in cell density, as more steps through the graph are required to move the same distance through a region of higher cell density.
-From the perspective of community detection algorithms, this effect "inflates" the high-density regions such that any internal substructure or noise is more likely to cause formation of subclusters.
-The resolution of clustering thus becomes dependent on the density of cells, which can occasionally be misleading if it overstates the heterogeneity in the data.
-
-```{r, echo=FALSE, eval=FALSE}
-set.seed(999)
-
-# For contributors who don't believe me, try this out.
-a <- matrix(rnorm(100000, 10), ncol=10)
-b <- matrix(rnorm(100, 0), ncol=10)
-x <- rbind(a, b)
-
-library(scran)
-g <- buildSNNGraph(x, transposed=TRUE, type="jaccard")
-out <- igraph::cluster_louvain(g)$membership
-table(out)
-```
-
-### Implementation
-
-To demonstrate, we use the `clusterCells()` function in `r Biocpkg("scran")` on PBMC dataset.
-All calculations are performed using the top PCs to take advantage of data compression and denoising.
-This function returns a vector containing cluster assignments for each cell in our `SingleCellExperiment` object.
-
-```{r}
-library(scran)
-nn.clusters <- clusterCells(sce.pbmc, use.dimred="PCA")
-table(nn.clusters)
-```
-
-We assign the cluster assignments back into our `SingleCellExperiment` object as a factor in the column metadata.
-This allows us to conveniently visualize the distribution of clusters in a $t$-SNE plot (Figure \@ref(fig:tsne-clust-graph)).
-
-```{r tsne-clust-graph, fig.cap="$t$-SNE plot of the 10X PBMC dataset, where each point represents a cell and is coloured according to the identity of the assigned cluster from graph-based clustering."}
-library(scater)
-colLabels(sce.pbmc) <- nn.clusters
-plotReducedDim(sce.pbmc, "TSNE", colour_by="label")
-```
-
-By default, `clusterCells()` uses the 10 nearest neighbors of each cell to construct a shared nearest neighbor graph.
-Two cells are connected by an edge if any of their nearest neighbors are shared,
-with the edge weight defined from the highest average rank of the shared neighbors [@xu2015identification].
-The Walktrap method from the `r CRANpkg("igraph")` package is then used to identify communities.
-
-If we wanted to explicitly specify all of these parameters, we would use the more verbose call below.
-This uses a `SNNGraphParam` object from the `r Biocpkg("bluster")` package to instruct `clusterCells()` to detect communities from a shared nearest-neighbor graph with the specified parameters.
-The appeal of this interface is that it allows us to easily switch to a different clustering algorithm by simply changing the `BLUSPARAM=` argument,
-as we will demonstrate later in the chapter.
-
-```{r}
-library(bluster)
-nn.clusters2 <- clusterCells(sce.pbmc, use.dimred="PCA",
- BLUSPARAM=SNNGraphParam(k=10, type="rank", cluster.fun="walktrap"))
-table(nn.clusters2)
-```
-
-We can also obtain the graph itself by specifying `full=TRUE` in the `clusterCells()` call.
-Doing so will return all intermediate structures that are used during clustering, including a graph object from the `r CRANpkg("igraph")` package.
-This graph can be visualized using a force-directed layout (Figure \@ref(fig:pbmc-force)), closely related to $t$-SNE and UMAP,
-though which of these is the most aesthetically pleasing is left to the eye of the beholder.
-
-```{r pbmc-force, message=FALSE, fig.cap="Force-directed layout for the shared nearest-neighbor graph of the PBMC dataset. Each point represents a cell and is coloured according to its assigned cluster identity."}
-nn.clust.info <- clusterCells(sce.pbmc, use.dimred="PCA", full=TRUE)
-nn.clust.info$objects$graph
-
-set.seed(11000)
-reducedDim(sce.pbmc, "force") <- igraph::layout_with_fr(nn.clust.info$objects$graph)
-plotReducedDim(sce.pbmc, colour_by="label", dimred="force")
-```
-
-In addition, the graph can be used to generate detailed diagnostics on the behavior of the graph-based clustering (`r link("using-graph-modularity", "OSCA.advanced")`).
-
-### Adjusting the parameters
-
-A graph-based clustering method has several key parameters:
-
-- How many neighbors are considered when constructing the graph.
-- What scheme is used to weight the edges.
-- Which community detection algorithm is used to define the clusters.
-
-One of the most important parameters is `k`, the number of nearest neighbors used to construct the graph.
-This controls the resolution of the clustering where higher `k` yields a more inter-connected graph and broader clusters.
-Users can exploit this by experimenting with different values of `k` to obtain a satisfactory resolution.
-
-```{r}
-# More resolved.
-clust.5 <- clusterCells(sce.pbmc, use.dimred="PCA", BLUSPARAM=NNGraphParam(k=5))
-table(clust.5)
-
-# Less resolved.
-clust.50 <- clusterCells(sce.pbmc, use.dimred="PCA", BLUSPARAM=NNGraphParam(k=50))
-table(clust.50)
-```
-
-```{r, echo=FALSE}
-stopifnot(length(unique(nn.clusters)) < length(unique(clust.5)))
-stopifnot(length(unique(nn.clusters)) > length(unique(clust.50)))
-```
-
-Further tweaking can be performed by changing the edge weighting scheme during graph construction.
-Setting `type="number"` will weight edges based on the number of nearest neighbors that are shared between two cells.
-Similarly, `type="jaccard"` will weight edges according to the Jaccard index of the two sets of neighbors.
-We can also disable weighting altogether by using a simple $k$-nearest neighbor graph, which is occasionally useful for downstream graph operations that do not support weights.
-
-```{r}
-clust.num <- clusterCells(sce.pbmc, use.dimred="PCA",
- BLUSPARAM=NNGraphParam(type="number"))
-table(clust.num)
-clust.jaccard <- clusterCells(sce.pbmc, use.dimred="PCA",
- BLUSPARAM=NNGraphParam(type="jaccard"))
-table(clust.jaccard)
-clust.none <- clusterCells(sce.pbmc, use.dimred="PCA",
- BLUSPARAM=KNNGraphParam())
-table(clust.none)
-```
-
-The community detection can be performed by using any of the algorithms provided by `r CRANpkg("igraph")`.
-We have already mentioned the Walktrap approach, but many others are available to choose from:
-
-```{r}
-clust.walktrap <- clusterCells(sce.pbmc, use.dimred="PCA",
- BLUSPARAM=NNGraphParam(cluster.fun="walktrap"))
-
-clust.louvain <- clusterCells(sce.pbmc, use.dimred="PCA",
- BLUSPARAM=NNGraphParam(cluster.fun="louvain"))
-
-clust.infomap <- clusterCells(sce.pbmc, use.dimred="PCA",
- BLUSPARAM=NNGraphParam(cluster.fun="infomap"))
-
-clust.fast <- clusterCells(sce.pbmc, use.dimred="PCA",
- BLUSPARAM=NNGraphParam(cluster.fun="fast_greedy"))
-
-clust.labprop <- clusterCells(sce.pbmc, use.dimred="PCA",
- BLUSPARAM=NNGraphParam(cluster.fun="label_prop"))
-
-clust.eigen <- clusterCells(sce.pbmc, use.dimred="PCA",
- BLUSPARAM=NNGraphParam(cluster.fun="leading_eigen"))
-```
-
-It is straightforward to compare two clustering strategies to see how they differ (`r link("comparing-different-clusterings", "OSCA.advanced")`).
-For example, Figure \@ref(fig:walktrap-v-others) suggests that Infomap yields finer clusters than Walktrap while fast-greedy yields coarser clusters.
-
-```{r walktrap-v-others, fig.asp=2, fig.cap="Number of cells assigned to combinations of cluster labels with different community detection algorithms in the PBMC dataset. Each entry of each heatmap represents a pair of labels, coloured proportionally to the log-number of cells with those labels."}
-library(pheatmap)
-
-# Using a large pseudo-count for a smoother color transition
-# between 0 and 1 cell in each 'tab'.
-tab <- table(paste("Infomap", clust.infomap),
- paste("Walktrap", clust.walktrap))
-ivw <- pheatmap(log10(tab+10), main="Infomap vs Walktrap",
- color=viridis::viridis(100), silent=TRUE)
-
-tab <- table(paste("Fast", clust.fast),
- paste("Walktrap", clust.walktrap))
-fvw <- pheatmap(log10(tab+10), main="Fast-greedy vs Walktrap",
- color=viridis::viridis(100), silent=TRUE)
-
-gridExtra::grid.arrange(ivw[[4]], fvw[[4]])
-```
-
-```{r, echo=FALSE}
-# Checking my comments above.
-stopifnot(length(unique(clust.fast)) < length(unique(clust.walktrap)))
-stopifnot(length(unique(clust.infomap)) > length(unique(clust.walktrap)))
-```
-
-Pipelines involving `r Biocpkg("scran")` default to rank-based weights followed by Walktrap clustering.
-In contrast, `r CRANpkg("Seurat")` uses Jaccard-based weights followed by Louvain clustering.
-Both of these strategies work well, and it is likely that the same could be said for many other combinations of weighting schemes and community detection algorithms.
-
-## Vector quantization with $k$-means
-
-### Background
-
-Vector quantization partitions observations into groups where each group is associated with a representative point, i.e., vector in the coordinate space.
-This is a type of clustering that primarily aims to compress data by replacing many points with a single representative.
-The representatives can then be treated as "samples" for further analysis, reducing the number of samples and computational work in later steps like, e.g., trajectory reconstruction [@ji2016tscan].
-This approach will also eliminate differences in cell density across the expression space, ensuring that the most abundant cell type does not dominate downstream results.
-
-$k$-means clustering is a classic vector quantization technique that divides cells into $k$ clusters.
-Each cell is assigned to the cluster with the closest centroid, which is done by minimizing the within-cluster sum of squares using a random starting configuration for the $k$ centroids.
-We usually set $k$ to a large value such as the square root of the number of cells to obtain fine-grained clusters.
-These are not meant to be interpreted directly, but rather, the centroids are used in downstream steps for faster computation.
-The main advantage of this approach lies in its speed, given the simplicity and ease of implementation of the algorithm.
-
-### Implementation
-
-We supply a `KmeansParam` object in `clusterCells()` to perform $k$-means clustering with the specified number of clusters in `centers=`.
-We again use our top PCs after setting the random seed to ensure that the results are reproducible.
-In general, the $k$-means clusters correspond to the visual clusters on the $t$-SNE plot in Figure \@ref(fig:tsne-clust-kmeans), though there are some divergences that are not observed in, say, Figure \@ref(fig:tsne-clust-graph).
-(This is at least partially due to the fact that $t$-SNE is itself graph-based and so will naturally agree more with a graph-based clustering strategy.)
-
-```{r tsne-clust-kmeans, fig.cap="$t$-SNE plot of the 10X PBMC dataset, where each point represents a cell and is coloured according to the identity of the assigned cluster from $k$-means clustering."}
-set.seed(100)
-clust.kmeans <- clusterCells(sce.pbmc, use.dimred="PCA",
- BLUSPARAM=KmeansParam(centers=10))
-table(clust.kmeans)
-
-colLabels(sce.pbmc) <- clust.kmeans
-plotReducedDim(sce.pbmc, "TSNE", colour_by="label")
-```
-
-If we were so inclined, we could obtain a "reasonable" choice of $k$ by computing the gap statistic using methods from the `r CRANpkg("cluster")` package.
-A more practical use of $k$-means is to deliberately set $k$ to a large value to achieve overclustering.
-This will forcibly partition cells inside broad clusters that do not have well-defined internal structure.
-For example, we might be interested in the change in expression from one "side" of a cluster to the other, but the lack of any clear separation within the cluster makes it difficult to separate with graph-based methods, even at the highest resolution.
-$k$-means has no such problems and will readily split these broad clusters for greater resolution.
-
-```{r tsne-clust-kmeans-best, fig.cap="$t$-SNE plot of the 10X PBMC dataset, where each point represents a cell and is coloured according to the identity of the assigned cluster from $k$-means clustering with $k=20$."}
-set.seed(100)
-clust.kmeans2 <- clusterCells(sce.pbmc, use.dimred="PCA",
- BLUSPARAM=KmeansParam(centers=20))
-table(clust.kmeans2)
-
-colLabels(sce.pbmc) <- clust.kmeans2
-plotTSNE(sce.pbmc, colour_by="label", text_by="label")
-```
-
-For larger datasets, we can use a variant of this approach named mini-batch $k$-means from the `r Biocpkg("mbkmeans")` package.
-At each iteration, we only update the cluster assignments and centroid positions for a small subset of the observations.
-This reduces memory usage and computational time - especially when not all of the observations are informative for convergence - and supports parallelization via `r Biocpkg("BiocParallel")`.
-Using this variant is as simple as switching to a `MbkmeansParam()` object in our `clusterCells()` call:
-
-```{r}
-set.seed(100)
-clust.mbkmeans <- clusterCells(sce.pbmc, use.dimred="PCA",
- BLUSPARAM=MbkmeansParam(centers=10))
-table(clust.mbkmeans)
-```
-
-### In two-step procedures
-
-By itself, $k$-means suffers from several shortcomings that reduce its appeal for obtaining interpretable clusters:
-
-- It implicitly favors spherical clusters of equal radius.
-This can lead to unintuitive partitionings on real datasets that contain groupings with irregular sizes and shapes.
-- The number of clusters $k$ must be specified beforehand and represents a hard cap on the resolution of the clustering..
-For example, setting $k$ to be below the number of cell types will always lead to co-clustering of two cell types, regardless of how well separated they are.
-In contrast, other methods like graph-based clustering will respect strong separation even if the relevant resolution parameter is set to a low value.
-- It is dependent on the randomly chosen initial coordinates.
-This requires multiple runs to verify that the clustering is stable.
-
-However, these concerns are less relevant when $k$-means is being used for vector quantization.
-In this application, $k$-means is used as a prelude to more sophisticated and interpretable - but computationally expensive - clustering algorithms.
-The `clusterCells()` function supports a "two-step" mode where $k$-means is initially used to obtain representative centroids that are subjected to graph-based clustering.
-Each cell is then placed in the same graph-based cluster that its $k$-means centroid was assigned to (Figure \@ref(fig:tsne-kmeans-graph-pbmc)).
-
-```{r tsne-kmeans-graph-pbmc, fig.cap="$t$-SNE plot of the PBMC dataset, where each point represents a cell and is coloured according to the identity of the assigned cluster from combined $k$-means/graph-based clustering."}
-# Setting the seed due to the randomness of k-means.
-set.seed(0101010)
-kgraph.clusters <- clusterCells(sce.pbmc, use.dimred="PCA",
- BLUSPARAM=TwoStepParam(
- first=KmeansParam(centers=1000),
- second=NNGraphParam(k=5)
- )
-)
-table(kgraph.clusters)
-
-plotTSNE(sce.pbmc, colour_by=I(kgraph.clusters))
-```
-
-The obvious benefit of this approach over direct graph-based clustering is the speed improvement.
-We avoid the need to identifying nearest neighbors for each cell and the construction of a large intermediate graph,
-while benefiting from the relative interpretability of graph-based clusters compared to those from $k$-means.
-This approach also mitigates the "inflation" effect discussed in Section \@ref(clustering-graph).
-Each centroid serves as a representative of a region of space that is roughly similar in volume,
-ameliorating differences in cell density that can cause (potentially undesirable) differences in resolution.
-
-The choice of the number of $k$-means clusters determines the trade-off between speed and fidelity.
-Larger values provide a more faithful representation of the underlying distribution of cells,
-at the cost of requiring more computational work by the second-step clustering procedure.
-Note that the second step operates on the centroids, so increasing `clusters=` may have further implications if the second-stage procedure is sensitive to the total number of input observations.
-For example, increasing the number of centroids would require an concomitant increase in `k=` (the number of neighbors in graph construction) to maintain the same level of resolution in the final output.
-
-## Hierarchical clustering
-
-### Background
-
-Hierarchical clustering is an old technique that arranges samples into a hierarchy based on their relative similarity to each other.
-Most implementations do so by joining the most similar samples into a new cluster, then joining similar clusters into larger clusters, and so on, until all samples belong to a single cluster.
-This process yields obtain a dendrogram that defines clusters with progressively increasing granularity.
-Variants of hierarchical clustering methods primarily differ in how they choose to perform the agglomerations.
-For example, complete linkage aims to merge clusters with the smallest maximum distance between their elements, while Ward's method aims to minimize the increase in within-cluster variance.
-
-In the context of scRNA-seq, the main advantage of hierarchical clustering lies in the production of the dendrogram.
-This is a rich summary that quantitatively captures the relationships between subpopulations at various resolutions.
-Cutting the dendrogram at high resolution is also guaranteed to yield clusters that are nested within those obtained at a low-resolution cut;
-this can be helpful for interpretation, as discussed in Section \@ref(subclustering).
-The dendrogram is also a natural representation of the data in situations where cells have descended from a relatively recent common ancestor.
-
-In practice, hierarchical clustering is too slow to be used for anything but the smallest scRNA-seq datasets.
-Most implementations require a cell-cell distance matrix that is prohibitively expensive to compute for a large number of cells.
-Greedy agglomeration is also likely to result in a quantitatively suboptimal partitioning (as defined by the agglomeration measure) at higher levels of the dendrogram when the number of cells and merge steps is high.
-Nonetheless, we will still demonstrate the application of hierarchical clustering here as it can be useful when combined with vector quantization techniques like $k$-means.
-
-### Implementation
-
-The PBMC dataset is too large to use directly in hierarchical clustering, requiring a two-step approach to compress the observations instead (Section \@ref(in-two-step-procedures)).
-For the sake of simplicity, we will demonstrate on the smaller 416B dataset instead.
-
-```{r, results='asis', echo=FALSE}
-extractFromPackage("lun-416b.Rmd", package="OSCA.workflows",
- chunk="dimensionality-reduction", objects="sce.416b")
-```
-
-```{r}
-sce.416b
-```
-
-We use a `HclustParam` object to instruct `clusterCells()` to perform hierarchical clustering on the top PCs.
-Specifically, it computes a cell-cell distance matrix using the top PCs and then applies Ward's minimum variance method to obtain a dendrogram.
-When visualized in Figure \@ref(fig:dend-416b), we see a clear split in the population caused by oncogene induction.
-While both Ward's method and the default complete linkage yield compact clusters, we prefer the former it is less affected by differences in variance between clusters.
-
-```{r dend-416b, fig.cap="Hierarchy of cells in the 416B data set after hierarchical clustering, where each leaf node is a cell that is coloured according to its oncogene induction status (red is induced, blue is control) and plate of origin (light or dark)."}
-hclust.416b <- clusterCells(sce.416b, use.dimred="PCA",
- BLUSPARAM=HclustParam(method="ward.D2"), full=TRUE)
-tree.416b <- hclust.416b$objects$hclust
-
-# Making a prettier dendrogram.
-library(dendextend)
-tree.416b$labels <- seq_along(tree.416b$labels)
-dend <- as.dendrogram(tree.416b, hang=0.1)
-
-combined.fac <- paste0(sce.416b$block, ".",
- sub(" .*", "", sce.416b$phenotype))
-labels_colors(dend) <- c(
- "20160113.wild"="blue",
- "20160113.induced"="red",
- "20160325.wild"="dodgerblue",
- "20160325.induced"="salmon"
-)[combined.fac][order.dendrogram(dend)]
-
-plot(dend)
-```
-
-```{r, echo=FALSE, results="hide", fig.show="hide"}
-set.seed(999)
-
-# For contributors who don't believe me, try this out.
-# Good luck trying to get a single cut that gets all three
-# clusters out with cutree!
-a1 <- matrix(rnorm(5000, 4, sd=2), ncol=10)
-a2 <- matrix(rnorm(5000, 2, sd=1.5), ncol=10)
-a3 <- matrix(rnorm(5000, 0, sd=1), ncol=10)
-x <- rbind(a1, a2, a3)
-
-dist.sim <- dist(x)
-tree.sim <- hclust(dist.sim, "ward.D2")
-tree.sim2 <- hclust(dist.sim)
-dend.sim <- as.dendrogram(tree.sim)
-dend.sim2 <- as.dendrogram(tree.sim2)
-
-cluster <- rep(1:3, each=500)
-labels_colors(dend.sim) <- c("blue", "red", "darkgreen")[cluster][order.dendrogram(dend.sim)]
-labels_colors(dend.sim2) <- c("blue", "red", "darkgreen")[cluster][order.dendrogram(dend.sim2)]
-
-par(mfrow=c(1,2))
-plot(dend.sim)
-plot(dend.sim2)
-
-table(cluster,
- dynamicTreeCut::cutreeDynamic(tree.sim, distM=as.matrix(dist.sim),
- minClusterSize=10))
-table(cluster,
- dynamicTreeCut::cutreeDynamic(tree.sim2, distM=as.matrix(dist.sim),
- minClusterSize=10))
-```
-
-To obtain explicit clusters, we "cut" the tree by removing internal branches such that every subtree represents a distinct cluster.
-This is most simply done by removing internal branches above a certain height of the tree, as performed by the `cutree()` function.
-A more sophisticated variant of this approach is implemented in the `r CRANpkg("dynamicTreeCut")` package,
-which uses the shape of the branches to obtain a better partitioning for complex dendrograms (Figure \@ref(fig:dend-cluster)).
-We enable this option by setting `cut.dynamic=TRUE`, with additional tweaking of the `deepSplit=` parameter to control the resolution of the resulting clusters.
-
-```{r dend-cluster, fig.cap="Hierarchy of cells in the 416B data set after hierarchical clustering, where each leaf node is a cell that is coloured according to its assigned cluster identity from a dynamic tree cut."}
-hclust.dyn <- clusterCells(sce.416b, use.dimred="PCA",
- BLUSPARAM=HclustParam(method="ward.D2", cut.dynamic=TRUE,
- cut.params=list(minClusterSize=10, deepSplit=1)))
-table(hclust.dyn)
-
-labels_colors(dend) <- as.integer(hclust.dyn)[order.dendrogram(dend)]
-plot(dend)
-```
-
-This generally corresponds well to the grouping of cells on a $t$-SNE plot (Figure \@ref(fig:tsne-416b)).
-Cluster 2 is split across two visual clusters in the plot but we attribute this to a distortion introduced by $t$-SNE,
-given that this cluster actually has the highest average silhouette width (`r link("silhouette-width", "OSCA.advanced")`).
-
-```{r tsne-416b, fig.cap="$t$-SNE plot of the 416B dataset, where each point represents a cell and is coloured according to the identity of the assigned cluster from hierarchical clustering."}
-colLabels(sce.416b) <- factor(hclust.dyn)
-plotReducedDim(sce.416b, "TSNE", colour_by="label")
-```
-
-```{r, echo=FALSE}
-set.seed(91019)
-fault <- 2L
-visual <- kmeans(reducedDim(sce.416b, "TSNE"), 4)
-xtab <- table(colLabels(sce.416b), visual$cluster)
-stopifnot(sum(xtab[fault,] > 0) == 2)
-stopifnot(sum(xtab[fault,]==0) == 2)
-
-sil <- cluster::silhouette(as.integer(hclust.dyn),
- dist = dist(reducedDim(sce.416b, "PCA")))
-sil.out <- vapply(split(sil[,3], sil[,1]), mean, 0)
-stopifnot(names(sil.out[which.max(sil.out)])=="2")
-```
-
-### In two-step procedures, again
-
-Returning to our PBMC example, we can use a two-step approach to perform hierarchical clustering on the representative centroids (Figure \@ref(fig:tsne-kmeans-hclust-pbmc)).
-This avoids the construction of a distance matrix across all cells for faster computation.
-
-```{r tsne-kmeans-hclust-pbmc, fig.cap="$t$-SNE plot of the PBMC dataset, where each point represents a cell and is coloured according to the identity of the assigned cluster from combined $k$-means/hierarchical clustering."}
-# Setting the seed due to the randomness of k-means.
-set.seed(1111)
-khclust.info <- clusterCells(sce.pbmc, use.dimred="PCA",
- BLUSPARAM=TwoStepParam(
- first=KmeansParam(centers=1000),
- second=HclustParam(method="ward.D2", cut.dynamic=TRUE,
- cut.param=list(deepSplit=3)) # for higher resolution.
- ),
- full=TRUE
-)
-table(khclust.info$clusters)
-
-plotTSNE(sce.pbmc, colour_by=I(khclust.info$clusters),
- text_by=I(khclust.info$clusters))
-```
-
-With a little bit of work, we can also examine the dendrogram constructed on the centroids (Figure \@ref(fig:dendroid-pbmc)).
-This provides a more quantitative visualization of the relative similarities between the different subpopulations.
-
-```{r dendroid-pbmc, fig.cap="Dendrogram of the $k$-mean centroids after hierarchical clustering in the PBMC dataset. Each leaf node represents a representative cluster of cells generated by $k$-mean clustering."}
-k.stats <- khclust.info$objects$first
-tree.pbmc <- khclust.info$objects$second$hclust
-
-m <- match(as.integer(tree.pbmc$labels), k.stats$cluster)
-final.clusters <- khclust.info$clusters[m]
-
-# TODO: expose scater color palette for easier re-use,
-# given that the default colors start getting recycled.
-dend <- as.dendrogram(tree.pbmc, hang=0.1)
-labels_colors(dend) <- as.integer(final.clusters)[order.dendrogram(dend)]
-
-plot(dend)
-```
-
-As an aside, the same approach can be used to speed up any clustering method based on a distance matrix.
-For example, we could subject our $k$-means centroids to clustering by affinity propagation [@frey2007clustering].
-In this procedure, each sample (i.e., centroid) chooses itself or another sample as its "exemplar",
-with the suitability of the choice dependent on the distance between the samples, other potential exemplars for each sample, and the other samples with the same chosen exemplar.
-Iterative updates of these choices yields a set of clusters where each cluster is defined from the samples assigned to the same exemplar (Figure \@ref(fig:tsne-kmeans-affinity-pbmc))
-Unlike hierarchical clustering, this does not provide a dendrogram, but it also avoids the extra complication of a tree cut -
-resolution is primarily controlled via the `q=` parameter, which defines the strength with which a sample considers itself as an exemplar and thus forms its own cluster.
-
-```{r, echo=FALSE}
-# Note technically necessary, but just throwing it in here
-# so that it gets picked up by the dependency scraper.
-library(apcluster)
-```
-
-```{r tsne-kmeans-affinity-pbmc, fig.cap="$t$-SNE plot of the PBMC dataset, where each point represents a cell and is coloured according to the identity of the assigned cluster from combined $k$-means/affinity propagation clustering."}
-# Setting the seed due to the randomness of k-means.
-set.seed(1111)
-kaclust.info <- clusterCells(sce.pbmc, use.dimred="PCA",
- BLUSPARAM=TwoStepParam(
- first=KmeansParam(centers=1000),
- second=AffinityParam(q=0.1) # larger q => more clusters
- ),
- full=TRUE
-)
-table(kaclust.info$clusters)
-
-plotTSNE(sce.pbmc, colour_by=I(kaclust.info$clusters),
- text_by=I(kaclust.info$clusters))
-```
-
-## Subclustering {#subclustering}
-
-Another simple approach to improving resolution is to repeat the feature selection and clustering _within_ a single cluster.
-This aims to select HVGs and PCs that are more relevant to internal structure, improving resolution by avoiding noise from unnecessary features.
-Subsetting also encourages clustering methods to separate cells according to more modest heterogeneity in the absence of distinct subpopulations.
-We demonstrate with a cluster of putative memory T cells from the PBMC dataset, identified according to several markers (Figure \@ref(fig:ccr7-dist-memory)).
-
-```{r ccr7-dist-memory, fig.cap="Distribution of log-normalized expression values for several T cell markers within each cluster in the 10X PBMC dataset. Each cluster is color-coded for convenience."}
-clust.full <- clusterCells(sce.pbmc, use.dimred="PCA")
-plotExpression(sce.pbmc, features=c("CD3E", "CCR7", "CD69", "CD44"),
- x=I(clust.full), colour_by=I(clust.full))
-
-# Repeating modelling and PCA on the subset.
-memory <- 6L
-sce.memory <- sce.pbmc[,clust.full==memory]
-dec.memory <- modelGeneVar(sce.memory)
-sce.memory <- denoisePCA(sce.memory, technical=dec.memory,
- subset.row=getTopHVGs(dec.memory, n=5000))
-```
-
-```{r, echo=FALSE}
-by.clust <- split(logcounts(sce.pbmc)["CCR7",], clust.full)
-by.clust.m <- vapply(by.clust, mean, 0)
-stopifnot(identical(unname(which.max(by.clust.m)), memory))
-```
-
-We apply graph-based clustering within this memory subset to obtain _CD4_^+^ and _CD8_^+^ subclusters (Figure \@ref(fig:cd48-memory-expression)).
-Admittedly, the expression of _CD4_ is so low that the change is rather modest, but the interpretation is clear enough.
-
-```{r cd48-memory-expression, fig.wide=TRUE, fig.asp=0.5, fig.cap="Distribution of _CD4_ and _CD8A_ log-normalized expression values within each cluster in the memory T cell subset of the 10X PBMC dataset."}
-g.memory <- buildSNNGraph(sce.memory, use.dimred="PCA")
-clust.memory <- igraph::cluster_walktrap(g.memory)$membership
-plotExpression(sce.memory, features=c("CD8A", "CD4"),
- x=I(factor(clust.memory)))
-```
-
-```{r, echo=FALSE}
-by.clust <- split(logcounts(sce.memory)["CD8A",], clust.memory)
-by.clust.cd8 <- vapply(by.clust, mean, 0)
-top.cd8 <- which.max(by.clust.cd8)
-stopifnot(length(by.clust.cd8) == 2L, all(by.clust.cd8[top.cd8] > by.clust.cd8[-top.cd8] + 0.8))
-
-by.clust <- split(logcounts(sce.memory)["CD4",], clust.memory)
-by.clust.cd4 <- vapply(by.clust, mean, 0)
-stopifnot(length(by.clust.cd4) == 2L, identical(which.min(by.clust.cd4), top.cd8))
-```
-
-For subclustering analyses, it is helpful to define a customized function that calls our desired algorithms to obtain a clustering from a given `SingleCellExperiment`.
-This function can then be applied multiple times on different subsets without having to repeatedly copy and modify the code for each subset.
-For example, `quickSubCluster()` loops over all subsets and executes this user-specified function to generate a list of `SingleCellExperiment` objects containing the subclustering results.
-(Of course, the downside is that this assumes that a similar analysis is appropriate for each subset.
-If different subsets require extensive reparametrization, copying the code may actually be more straightforward.)
-
-```{r}
-set.seed(1000010)
-subcluster.out <- quickSubCluster(sce.pbmc, groups=clust.full,
- prepFUN=function(x) { # Preparing the subsetted SCE for clustering.
- dec <- modelGeneVar(x)
- input <- denoisePCA(x, technical=dec,
- subset.row=getTopHVGs(dec, prop=0.1),
- BSPARAM=BiocSingular::IrlbaParam())
- },
- clusterFUN=function(x) { # Performing the subclustering in the subset.
- g <- buildSNNGraph(x, use.dimred="PCA", k=20)
- igraph::cluster_walktrap(g)$membership
- }
-)
-
-# One SingleCellExperiment object per parent cluster:
-names(subcluster.out)
-
-# Looking at the subclustering for one example:
-table(subcluster.out[[1]]$subcluster)
-```
-
-Subclustering is a general and conceptually straightforward procedure for increasing resolution.
-It can also simplify the interpretation of the subclusters, which only need to be considered in the context of the parent cluster's identity - for example, we did not have to re-identify the cells in cluster `r memory` as T cells.
-However, this is a double-edged sword as it is difficult for practitioners to consider the uncertainty of identification for parent clusters when working with deep nesting.
-If cell types or states span cluster boundaries, conditioning on the putative cell type identity of the parent cluster can encourage the construction of a "house of cards" of cell type assignments, e.g., where a subcluster of one parent cluster is actually contamination from a cell type in a separate parent cluster.
-
-## Session Info {-}
-
-```{r sessionInfo, echo=FALSE, results='asis'}
-prettySessionInfo()
-```
diff --git a/inst/book/clustering.html b/inst/book/clustering.html
new file mode 100644
index 0000000..830cee9
--- /dev/null
+++ b/inst/book/clustering.html
@@ -0,0 +1,9 @@
+
+
+
+
+
+
+
+
+
diff --git a/inst/book/dimensionality-reduction.html b/inst/book/dimensionality-reduction.html
new file mode 100644
index 0000000..fde2ca1
--- /dev/null
+++ b/inst/book/dimensionality-reduction.html
@@ -0,0 +1,9 @@
+
+
+
+
+
+
+
+
+
diff --git a/inst/book/feature-selection.Rmd b/inst/book/feature-selection.Rmd
deleted file mode 100644
index 8a28748..0000000
--- a/inst/book/feature-selection.Rmd
+++ /dev/null
@@ -1,275 +0,0 @@
----
-output: html_document
-bibliography: ref.bib
----
-
-# Feature selection
-
-```{r setup, echo=FALSE, results="asis"}
-library(rebook)
-chapterPreamble()
-```
-
-## Motivation
-
-We often use scRNA-seq data in exploratory analyses to characterize heterogeneity across cells.
-Procedures like clustering and dimensionality reduction compare cells based on their gene expression profiles, which involves aggregating per-gene differences into a single (dis)similarity metric between a pair of cells.
-The choice of genes to use in this calculation has a major impact on the behavior of the metric and the performance of downstream methods.
-We want to select genes that contain useful information about the biology of the system while removing genes that contain random noise.
-This aims to preserve interesting biological structure without the variance that obscures that structure, and to reduce the size of the data to improve computational efficiency of later steps.
-
-The simplest approach to feature selection is to select the most variable genes based on their expression across the population.
-This assumes that genuine biological differences will manifest as increased variation in the affected genes, compared to other genes that are only affected by technical noise or a baseline level of "uninteresting" biological variation (e.g., from transcriptional bursting).
-Several methods are available to quantify the variation per gene and to select an appropriate set of highly variable genes (HVGs).
-We will discuss these below using the 10X PBMC dataset for demonstration:
-
-```{r, echo=FALSE, results="asis"}
-extractFromPackage("tenx-unfiltered-pbmc4k.Rmd", package="OSCA.workflows",
- chunk="normalization", objects="sce.pbmc")
-```
-
-```{r}
-sce.pbmc
-```
-
-As well as the 416B dataset:
-
-```{r, echo=FALSE, results="asis"}
-extractFromPackage("lun-416b.Rmd", package="OSCA.workflows",
- chunk="normalization", objects="sce.416b")
-```
-
-```{r}
-sce.416b
-```
-
-## Quantifying per-gene variation
-
-The simplest approach to quantifying per-gene variation is to compute the variance of the log-normalized expression values (i.e., "log-counts" ) for each gene across all cells [@lun2016step].
-The advantage of this approach is that the feature selection is based on the same log-values that are used for later downstream steps.
-In particular, genes with the largest variances in log-values will contribute most to the Euclidean distances between cells during procedures like clustering and dimensionality reduction.
-By using log-values here, we ensure that our quantitative definition of heterogeneity is consistent throughout the entire analysis.
-
-Calculation of the per-gene variance is simple but feature selection requires modelling of the mean-variance relationship.
-The log-transformation is not a variance stabilizing transformation in most cases,
-which means that the total variance of a gene is driven more by its abundance than its underlying biological heterogeneity.
-To account for this effect, we use the `modelGeneVar()` function to fit a trend to the variance with respect to abundance across all genes (Figure \@ref(fig:trend-plot-pbmc)).
-
-```{r trend-plot-pbmc, fig.cap="Variance in the PBMC data set as a function of the mean. Each point represents a gene while the blue line represents the trend fitted to all genes."}
-library(scran)
-dec.pbmc <- modelGeneVar(sce.pbmc)
-
-# Visualizing the fit:
-fit.pbmc <- metadata(dec.pbmc)
-plot(fit.pbmc$mean, fit.pbmc$var, xlab="Mean of log-expression",
- ylab="Variance of log-expression")
-curve(fit.pbmc$trend(x), col="dodgerblue", add=TRUE, lwd=2)
-```
-
-At any given abundance, we assume that the variation in expression for most genes is driven by uninteresting processes like sampling noise.
-Under this assumption, the fitted value of the trend at any given gene's abundance represents an estimate of its uninteresting variation, which we call the technical component.
-We then define the biological component for each gene as the difference between its total variance and the technical component.
-This biological component represents the "interesting" variation for each gene and can be used as the metric for HVG selection.
-
-```{r}
-# Ordering by most interesting genes for inspection.
-dec.pbmc[order(dec.pbmc$bio, decreasing=TRUE),]
-```
-
-(Careful readers will notice that some genes have negative biological components, which have no obvious interpretation and can be ignored in most applications.
-They are inevitable when fitting a trend to the per-gene variances as approximately half of the genes will lie below the trend.)
-
-Strictly speaking, the interpretation of the fitted trend as the technical component assumes that the expression profiles of most genes are dominated by random technical noise.
-In practice, all expressed genes will exhibit some non-zero level of biological variability due to events like transcriptional bursting.
-Thus, it would be more appropriate to consider these estimates as technical noise plus "uninteresting" biological variation,
-under the assumption that most genes do not participate in the processes driving interesting heterogeneity across the population.
-
-## Quantifying technical noise {#sec:spikeins}
-
-The assumption in Section \@ref(quantifying-per-gene-variation) may be problematic in rare scenarios where many genes at a particular abundance are affected by a biological process.
-For example, strong upregulation of cell type-specific genes may result in an enrichment of HVGs at high abundances.
-This would inflate the fitted trend in that abundance interval and compromise the detection of the relevant genes.
-We can avoid this problem by fitting a mean-dependent trend to the variance of the spike-in transcripts (Figure \@ref(fig:spike-416b)), if they are available.
-The premise here is that spike-ins should not be affected by biological variation, so the fitted value of the spike-in trend should represent a better estimate of the technical component for each gene.
-
-```{r spike-416b, fig.cap="Variance in the 416B data set as a function of the mean. Each point represents a gene (black) or spike-in transcript (red) and the blue line represents the trend fitted to all spike-ins."}
-dec.spike.416b <- modelGeneVarWithSpikes(sce.416b, "ERCC")
-dec.spike.416b[order(dec.spike.416b$bio, decreasing=TRUE),]
-
-plot(dec.spike.416b$mean, dec.spike.416b$total, xlab="Mean of log-expression",
- ylab="Variance of log-expression")
-fit.spike.416b <- metadata(dec.spike.416b)
-points(fit.spike.416b$mean, fit.spike.416b$var, col="red", pch=16)
-curve(fit.spike.416b$trend(x), col="dodgerblue", add=TRUE, lwd=2)
-```
-
-In the absence of spike-in data, one can attempt to create a trend by making some distributional assumptions about the noise.
-For example, UMI counts typically exhibit near-Poisson variation if we only consider technical noise from library preparation and sequencing.
-This can be used to construct a mean-variance trend in the log-counts (Figure \@ref(fig:tech-pbmc)) with the `modelGeneVarByPoisson()` function.
-Note the increased residuals of the high-abundance genes, which can be interpreted as the amount of biological variation that was assumed to be "uninteresting" when fitting the gene-based trend in Figure \@ref(fig:trend-plot-pbmc).
-
-```{r tech-pbmc, fig.cap="Variance of normalized log-expression values for each gene in the PBMC dataset, plotted against the mean log-expression. The blue line represents represents the mean-variance relationship corresponding to Poisson noise."}
-set.seed(0010101)
-dec.pois.pbmc <- modelGeneVarByPoisson(sce.pbmc)
-dec.pois.pbmc <- dec.pois.pbmc[order(dec.pois.pbmc$bio, decreasing=TRUE),]
-head(dec.pois.pbmc)
-
-plot(dec.pois.pbmc$mean, dec.pois.pbmc$total, pch=16, xlab="Mean of log-expression",
- ylab="Variance of log-expression")
-curve(metadata(dec.pois.pbmc)$trend(x), col="dodgerblue", add=TRUE)
-```
-
-Interestingly, trends based purely on technical noise tend to yield large biological components for highly-expressed genes.
-This often includes so-called "house-keeping" genes coding for essential cellular components such as ribosomal proteins, which are considered uninteresting for characterizing cellular heterogeneity.
-These observations suggest that a more accurate noise model does not necessarily yield a better ranking of HVGs, though one should keep an open mind - house-keeping genes are regularly DE in a variety of conditions [@glare2002betaactin;@nazari2015gapdh;@guimaraes2016patterns], and the fact that they have large biological components indicates that there is strong variation across cells that may not be completely irrelevant.
-
-## Handling batch effects {#variance-batch}
-
-Data containing multiple batches will often exhibit batch effects - see `r link("integrating-datasets","OSCA.multisample")` for more details.
-We are usually not interested in HVGs that are driven by batch effects; instead, we want to focus on genes that are highly variable within each batch.
-This is naturally achieved by performing trend fitting and variance decomposition separately for each batch.
-We demonstrate this approach by treating each plate (`block`) in the 416B dataset as a different batch, using the `modelGeneVarWithSpikes()` function.
-(The same argument is available in all other variance-modelling functions.)
-
-```{r}
-dec.block.416b <- modelGeneVarWithSpikes(sce.416b, "ERCC", block=sce.416b$block)
-head(dec.block.416b[order(dec.block.416b$bio, decreasing=TRUE),1:6])
-```
-
-The use of a batch-specific trend fit is useful as it accommodates differences in the mean-variance trends between batches.
-This is especially important if batches exhibit systematic technical differences, e.g., differences in coverage or in the amount of spike-in RNA added.
-In this case, there are only minor differences between the trends in Figure \@ref(fig:blocked-fit), which indicates that the experiment was tightly replicated across plates.
-The analysis of each plate yields estimates of the biological and technical components for each gene, which are averaged across plates to take advantage of information from multiple batches.
-
-```{r blocked-fit, fig.asp=0.5, fig.width=10, fig.cap="Variance in the 416B data set as a function of the mean after blocking on the plate of origin. Each plot represents the results for a single plate, each point represents a gene (black) or spike-in transcript (red) and the blue line represents the trend fitted to all spike-ins."}
-par(mfrow=c(1,2))
-blocked.stats <- dec.block.416b$per.block
-for (i in colnames(blocked.stats)) {
- current <- blocked.stats[[i]]
- plot(current$mean, current$total, main=i, pch=16, cex=0.5,
- xlab="Mean of log-expression", ylab="Variance of log-expression")
- curfit <- metadata(current)
- points(curfit$mean, curfit$var, col="red", pch=16)
- curve(curfit$trend(x), col='dodgerblue', add=TRUE, lwd=2)
-}
-```
-
-Alternatively, we might consider using a linear model to account for batch effects and other unwanted factors of variation.
-This is more flexible as it can handle multiple factors and continuous covariates, though it is less accurate than `block=` in the special case of a multi-batch design.
-See `r link("handling-covariates-with-linear-models", "OSCA.advanced")` for more details.
-
-As an aside, the wave-like shape observed above is typical of the mean-variance trend for log-expression values.
-(The same wave is present but much less pronounced for UMI data.)
-A linear increase in the variance is observed as the mean increases from zero, as larger variances are obviously possible when the counts are not all equal to zero.
-In contrast, the relative contribution of sampling noise decreases at high abundances, resulting in a downward trend.
-The peak represents the point at which these two competing effects cancel each other out.
-
-## Selecting highly variable genes {#hvg-selection}
-
-Once we have quantified the per-gene variation, the next step is to select the subset of HVGs to use in downstream analyses.
-A larger subset will reduce the risk of discarding interesting biological signal by retaining more potentially relevant genes, at the cost of increasing noise from irrelevant genes that might obscure said signal.
-It is difficult to determine the optimal trade-off for any given application as noise in one context may be useful signal in another.
-For example, heterogeneity in T cell activation responses is an interesting phenomena [@richard2018tcell] but may be irrelevant noise in studies that only care about distinguishing the major immunophenotypes.
-
-The most obvious selection strategy is to take the top $n$ genes with the largest values for the relevant variance metric.
-The main advantage of this approach is that the user can directly control the number of genes retained, which ensures that the computational complexity of downstream calculations is easily predicted.
-For `modelGeneVar()` and `modelGeneVarWithSpikes()`, we would select the genes with the largest biological components.
-This is conveniently done for us via `getTopHVgs()`, as shown below with $n=1000$.
-
-```{r}
-# Taking the top 1000 genes here:
-hvg.pbmc.var <- getTopHVGs(dec.pbmc, n=1000)
-str(hvg.pbmc.var)
-```
-
-The choice of $n$ also has a fairly straightforward biological interpretation.
-Recall our trend-fitting assumption that most genes do not exhibit biological heterogeneity; this implies that they are not differentially expressed between cell types or states in our population.
-If we quantify this assumption into a statement that, e.g., no more than 5% of genes are differentially expressed, we can naturally set $n$ to 5% of the number of genes.
-In practice, we usually do not know the proportion of DE genes beforehand so this interpretation just exchanges one unknown for another.
-Nonetheless, it is still useful as it implies that we should lower $n$ for less heterogeneous datasets, retaining most of the biological signal without unnecessary noise from irrelevant genes.
-Conversely, more heterogeneous datasets should use larger values of $n$ to preserve secondary factors of variation beyond those driving the most obvious HVGs.
-
-The main disadvantage of this approach that it turns HVG selection into a competition between genes, whereby a subset of very highly variable genes can push other informative genes out of the top set.
-This can be problematic for analyses of highly heterogeneous populations if the loss of important markers prevents the resolution of certain subpopulations.
-In the most extreme example, consider a situation where a single subpopulation is very different from the others.
-In such cases, the top set will be dominated by differentially expressed genes involving that distinct subpopulation, compromising resolution of heterogeneity between the other populations.
-(This can be recovered with a nested analysis, as discussed in Section \@ref(subclustering), but we would prefer to avoid the problem in the first place.)
-
-
-
-Another potential concern with this approach is the fact that the choice of $n$ is fairly arbitrary, with any value from 500 to 5000 considered "reasonable".
-We have chosen $n=1000$ in the code above though there is no particular _a priori_ reason for doing so.
-Our recommendation is to simply pick an arbitrary $n$ and proceed with the rest of the analysis, with the intention of testing other choices later, rather than spending much time worrying about obtaining the "optimal" value.
-Alternatively, we may pick one of the other selection strategies discussed in `r link("more-hvg-selection-strategies", "OSCA.advanced")`.
-
-## Putting it all together {#feature-selection-subsetting}
-
-The code chunk below will select the top 10% of genes with the highest biological components.
-
-```{r}
-dec.pbmc <- modelGeneVar(sce.pbmc)
-chosen <- getTopHVGs(dec.pbmc, prop=0.1)
-str(chosen)
-```
-
-We then have several options to enforce our HVG selection on the rest of the analysis.
-
-- We can subset the `SingleCellExperiment` to only retain our selection of HVGs.
-This ensures that downstream methods will only use these genes for their calculations.
-The downside is that the non-HVGs are discarded from the new `SingleCellExperiment`, making it slightly more inconvenient to interrogate the full dataset for interesting genes that are not HVGs.
-
- ```{r}
- sce.pbmc.hvg <- sce.pbmc[chosen,]
- dim(sce.pbmc.hvg)
- ```
-
-- We can keep the original `SingleCellExperiment` object and specify the genes to use for downstream functions via an extra argument like `subset.row=`.
-This is useful if the analysis uses multiple sets of HVGs at different steps, whereby one set of HVGs can be easily swapped for another in specific steps.
-
- ```{r}
- # Performing PCA only on the chosen HVGs.
- library(scater)
- sce.pbmc <- runPCA(sce.pbmc, subset_row=chosen)
- reducedDimNames(sce.pbmc)
- ```
-
- This approach is facilitated by the `rowSubset()` utility,
- which allows us to easily store one or more sets of interest in our `SingleCellExperiment`.
- By doing so, we avoid the need to keep track of a separate `chosen` variable
- and ensure that our HVG set is synchronized with any downstream row subsetting of `sce.pbmc`.
-
- ```{r}
- rowSubset(sce.pbmc) <- chosen # stored in the default 'subset'.
- rowSubset(sce.pbmc, "HVGs.more") <- getTopHVGs(dec.pbmc, prop=0.2)
- rowSubset(sce.pbmc, "HVGs.less") <- getTopHVGs(dec.pbmc, prop=0.3)
- colnames(rowData(sce.pbmc))
- ```
-
- It can be inconvenient to repeatedly specify the desired feature set across steps,
- so some downstream functions will automatically subset to the default `rowSubset()` if present in the `SingleCellExperiment`.
- However, we find that it is generally safest to be explicit about which set is being used for a particular step.
-
-- We can have our cake and eat it too by (ab)using the "alternative Experiment" system in the `SingleCellExperiment` class.
-Initially designed for storing alternative features like spike-ins or antibody tags, we can instead use it to hold our full dataset while we perform our downstream operations conveniently on the HVG subset.
-This avoids book-keeping problems in long analyses when the original dataset is not synchronized with the HVG subsetted data.
-
- ```{r}
- # Recycling the class above.
- altExp(sce.pbmc.hvg, "original") <- sce.pbmc
- altExpNames(sce.pbmc.hvg)
-
- # No need for explicit subset_row= specification in downstream operations.
- sce.pbmc.hvg <- runPCA(sce.pbmc.hvg)
-
- # Recover original data:
- sce.pbmc.original <- altExp(sce.pbmc.hvg, "original", withColData=TRUE)
- ```
-
-## Session Info {-}
-
-```{r sessionInfo, echo=FALSE, results='asis'}
-prettySessionInfo()
-```
diff --git a/inst/book/feature-selection.html b/inst/book/feature-selection.html
new file mode 100644
index 0000000..ff5edbb
--- /dev/null
+++ b/inst/book/feature-selection.html
@@ -0,0 +1,10 @@
+
+
+
+
+
+
+
+
+
+
diff --git a/inst/book/index.Rmd b/inst/book/index.Rmd
deleted file mode 100644
index 3224e6a..0000000
--- a/inst/book/index.Rmd
+++ /dev/null
@@ -1,31 +0,0 @@
----
-title: "Basics of Single-Cell Analysis with Bioconductor"
-documentclass: book
-bibliography: ref.bib
-biblio-style: apalike
-link-citations: yes
-favicon: "`r rebook::BiocFavicon()`"
-github-repo: "OSCA-source/OSCA.basic"
----
-
-```{r, echo=FALSE}
-.gh_url <- file.path('https://github.com', rmarkdown::metadata[['github-repo']])
-```
-
----
-date: "`r rebook::openingDetails(Source=.gh_url)`"
-url: "`r .gh_url`"
----
-
-# Welcome {-}
-
-
-
-This site contains the basic analysis chapters for the `r BiocStyle::Biocbook("OSCA", label='__"Orchestrating Single-Cell Analysis with Bioconductor"__ book')`.
-```{r, results="asis", echo=FALSE}
-text <- packageDescription("OSCA.basic")$Description
-remainders <- strsplit(text, "\\.\\s")[[1]][-1]
-remainders <- sub("^\\s+", "", remainders)
-remainders <- paste(remainders, collapse=". ")
-cat(remainders, sep="\n")
-```
diff --git a/inst/book/index.html b/inst/book/index.html
new file mode 100644
index 0000000..185bdb7
--- /dev/null
+++ b/inst/book/index.html
@@ -0,0 +1,9 @@
+
+
+
+
+
+
+
+
+
diff --git a/inst/book/marker-detection.Rmd b/inst/book/marker-detection.Rmd
deleted file mode 100644
index 6a0f5ed..0000000
--- a/inst/book/marker-detection.Rmd
+++ /dev/null
@@ -1,348 +0,0 @@
----
-output:
- html_document
-bibliography: ref.bib
----
-
-# Marker gene detection {#marker-detection}
-
-```{r setup, echo=FALSE, results="asis"}
-library(rebook)
-chapterPreamble()
-```
-
-## Motivation
-
-To interpret our clustering results from Chapter \@ref(clustering), we identify the genes that drive separation between clusters.
-These marker genes allow us to assign biological meaning to each cluster based on their functional annotation.
-In the simplest case, we have _a priori_ knowledge of the marker genes associated with particular cell types, allowing us to treat the clustering as a proxy for cell type identity.
-The same principle can be applied to discover more subtle differences between clusters (e.g., changes in activation or differentiation state) based on the behavior of genes in the affected pathways.
-
-The most straightforward approach to marker gene detection involves testing for differential expression between clusters.
-If a gene is strongly DE between clusters, it is likely to have driven the separation of cells in the clustering algorithm.
-Several methods are available to quantify the differences in expression profiles between clusters and obtain a single ranking of genes for each cluster.
-We will demonstrate some of these choices in this chapter using the 10X PBMC dataset:
-
-```{r, results='asis', echo=FALSE}
-extractFromPackage("tenx-unfiltered-pbmc4k.Rmd", package="OSCA.workflows",
- chunk="clustering", objects="sce.pbmc")
-```
-
-```{r}
-sce.pbmc
-```
-
-## Scoring markers by pairwise comparisons
-
-```{r, echo=FALSE}
-target <- "5"
-```
-
-Our general strategy is to compare each pair of clusters and compute scores quantifying the differences in the expression distributions between clusters.
-The scores for all pairwise comparisons involving a particular cluster are then consolidated into a single `DataFrame` for that cluster.
-The `scoreMarkers()` function from `r Biocpkg("scran")` returns a list of `DataFrame`s where each `DataFrame` corresponds to a cluster and each row of the `DataFrame` corresponds to a gene.
-In the `DataFrame` for cluster $X$, the columns contain the `self.average`, the mean log-expression in $X$;
-`other.average`, the grand mean across all other clusters;
-`self.detected`, the proportion of cells with detected expression in $X$;
-`other.detected`, the mean detected proportion across all other clusters;
-and finally, a variety of effect size summaries generated from all pairwise comparisons involving $X$.
-
-```{r}
-library(scran)
-marker.info <- scoreMarkers(sce.pbmc, colLabels(sce.pbmc))
-marker.info
-colnames(marker.info[["5"]]) # statistics for cluster 5.
-```
-
-For each cluster, we can then rank candidate markers based on one of these effect size summaries.
-We demonstrate below with the mean AUC for cluster `r target`, which probably contains NK cells based on the top genes in Figure \@ref(fig:pbmc-scored-markers-5) (and no _CD3E_ expression).
-The next section will go into more detail on the differences between the various columns.
-
-```{r pbmc-scored-markers-5, fig.height=10, fig.width=8, fig.cap=sprintf("Distribution of expression values across clusters for the top potential marker genes (as determined by the mean AUC) for cluster %s in the PBMC dataset.", target)}
-chosen <- marker.info[["5"]]
-ordered <- chosen[order(chosen$mean.AUC, decreasing=TRUE),]
-head(ordered[,1:4]) # showing basic stats only, for brevity.
-
-library(scater)
-plotExpression(sce.pbmc, features=head(rownames(ordered)),
- x="label", colour_by="label")
-```
-
-```{r, echo=FALSE}
-# Checking the code and text are synced.
-stopifnot(identical(chosen[1,], marker.info[[target]][1,]))
-
-# Also noting that CD3 is not highly expressed.
-stopifnot(which("CD3E"==rownames(ordered)) > 10000)
-stopifnot(all(c("GNLY", "NKG7") %in% head(rownames(ordered))))
-```
-
-We deliberately use pairwise comparisons rather than comparing each cluster to the average of all other cells.
-The latter approach is sensitive to the population composition, which introduces an element of unpredictability to the marker sets due to variation in cell type abundances.
-(In the worst case, the presence of one subpopulation containing a majority of the cells will drive the selection of top markers for every other cluster, pushing out useful genes that can distinguish between the smaller subpopulations.)
-Moreover, pairwise comparisons naturally provide more information to interpret of the utility of a marker, e.g., by providing log-fold changes to indicate which clusters are distinguished by each gene (Section \@ref(obtaining-the-full-effects)).
-
-Previous editions of this chapter used $p$-values from the tests corresponding to each effect size, e.g., Welch's $t$-test, the Wilcoxon ranked sum test.
-While this is fine for ranking genes, the $p$-values themselves are statistically flawed and are of little use for inference -
-see `r link("p-value-invalidity", "OSCA.advanced")` for more details.
-The `scoreMarkers()` function simplifies the marker detection procedure by omitting the $p$-values altogether, instead focusing on the underlying effect sizes.
-
-## Effect sizes for pairwise comparisons
-
-In the context of marker detection, the area under the curve (AUC) quantifies our ability to distinguish between two distributions in a pairwise comparison.
-The AUC represents the probability that a randomly chosen observation from our cluster of interest is greater than a randomly chosen observation from the other cluster.
-A value of 1 corresponds to upregulation, where all values of our cluster of interest are greater than any value from the other cluster;
-a value of 0.5 means that there is no net difference in the location of the distributions;
-and a value of 0 corresponds to downregulation.
-The AUC is closely related to the $U$ statistic in the Wilcoxon ranked sum test (a.k.a., Mann-Whitney U-test).
-
-```{r}
-auc.only <- chosen[,grepl("AUC", colnames(chosen))]
-auc.only[order(auc.only$mean.AUC,decreasing=TRUE),]
-```
-
-Cohen's $d$ is a standardized log-fold change where the difference in the mean log-expression between groups is scaled by the average standard deviation across groups.
-In other words, it is the number of standard deviations that separate the means of the two groups.
-The interpretation is similar to the log-fold change; positive values indicate that the gene is upregulated in our cluster of interest,
-negative values indicate downregulation and values close to zero indicate that there is little difference.
-Cohen's $d$ is roughly analogous to the $t$-statistic in various two-sample $t$-tests.
-
-```{r}
-cohen.only <- chosen[,grepl("logFC.cohen", colnames(chosen))]
-cohen.only[order(cohen.only$mean.logFC.cohen,decreasing=TRUE),]
-```
-
-Finally, we also compute the log-fold change in the proportion of cells with detected expression between clusters.
-This ignores any information about the magnitude of expression, only considering whether any expression is detected at all.
-Again, positive values indicate that a greater proportion of cells express the gene in our cluster of interest compared to the other cluster.
-Note that a pseudo-count is added to avoid undefined log-fold changes when no cells express the gene in either group.
-
-```{r}
-detect.only <- chosen[,grepl("logFC.detected", colnames(chosen))]
-detect.only[order(detect.only$mean.logFC.detected,decreasing=TRUE),]
-```
-
-The AUC or Cohen's $d$ is usually the best choice for general purpose marker detection, as they are effective regardless of the magnitude of the expression values.
-The log-fold change in the detected proportion is specifically useful for identifying binary changes in expression.
-See `r link("properties-of-each-effect-size", "OSCA.advanced")` for more information about the practical differences between the effect sizes.
-
-## Summarizing pairwise effects
-
-In a dataset with $N$ clusters, each cluster is associated with $N-1$ values for each type of effect size described in the previous section.
-To simplify interpretation, we summarize the effects for each cluster into some key statistics such as the mean and median.
-Each summary statistic has a different interpretation when used for ranking:
-
-- The most obvious summary statistic is the mean.
-For cluster $X$, a large mean effect size (>0 for the log-fold changes, >0.5 for the AUCs) indicates that the gene is upregulated in $X$ compared to the average of the other groups.
-- Another summary statistic is the median, where a large value indicates that the gene is upregulated in $X$ compared to most (>50%) other clusters.
-The median provides greater robustness to outliers than the mean, which may or may not be desirable.
-On one hand, the median avoids an inflated effect size if only a minority of comparisons have large effects;
-on the other hand, it will also overstate the effect size by ignoring a minority of comparisons that have opposing effects.
-- The minimum value (`min.*`) is the most stringent summary for identifying upregulated genes, as a large value indicates that the gene is upregulated in $X$ compared to _all_ other clusters.
-Conversely, if the minimum is small (<0 for the log-fold changes, <0.5 for the AUCs), we can conclude that the gene is downregulated in $X$ compared to at least one other cluster.
-- The maximum value (`max.*`) is the least stringent summary for identifying upregulated genes, as a large value can be obtained if there is strong upregulation in $X$ compared to _any_ other cluster.
-Conversely, if the maximum is small, we can conclude that the gene is downregulated in $X$ compared to all other clusters.
-- The minimum rank, a.k.a., "min-rank" (`rank.*`) is the smallest rank of each gene across all pairwise comparisons.
-Specifically, genes are ranked _within_ each pairwise comparison based on decreasing effect size, and then the smallest rank _across_ all comparisons is reported for each gene.
-If a gene has a small min-rank, we can conclude that it is one of the top upregulated genes in at least one comparison of $X$ to another cluster.
-
-```{r, echo=FALSE}
-target <- "12"
-```
-
-Each of these summaries is computed for each effect size, for each gene, and for each cluster.
-Our next step is to choose one of these summary statistics for one of the effect sizes and to use it to rank the rows of the `DataFrame`.
-The choice of summary determines the stringency of the marker selection strategy, i.e., how many other clusters must we differ from?
-For identifying upregulated genes, ranking by the minimum is the most stringent and the maximum is the least stringent;
-the mean and median fall somewhere in between and are reasonable defaults for most applications.
-The example below uses the median Cohen's $d$ to obtain a ranking of upregulated markers for cluster `r target` (Figure \@ref(fig:pbmc-scored-markers-again)), which probably contains monocytes.
-
-```{r pbmc-scored-markers-again, fig.height=10, fig.width=8, fig.cap=sprintf("Distribution of expression values across clusters for the top potential marker genes (as determined by the median Cohen's $d$) for cluster %s in the PBMC dataset.", target)}
-chosen <- marker.info[["12"]] # using another cluster, for some variety.
-ordered <- chosen[order(chosen$median.logFC.cohen,decreasing=TRUE),]
-head(ordered[,1:4]) # showing basic stats only, for brevity.
-plotExpression(sce.pbmc, features=head(rownames(ordered)),
- x="label", colour_by="label")
-```
-
-```{r, echo=FALSE}
-# Checking the code and text are synced.
-stopifnot(identical(chosen[1,], marker.info[[target]][1,]))
-
-# Higher CD14 than CD16.
-stopifnot(which("CD14"==rownames(ordered)) < 50)
-stopifnot(which("FCGR3A"==rownames(ordered)) > 500)
-stopifnot(all(c("LYZ", "S100A9", "S100A8") %in% head(rownames(ordered))))
-```
-
-On some occasions, ranking by the minimum can be highly effective as it yields a concise set of highly cluster-specific markers.
-However, any gene that is expressed at the same level in two or more clusters will simply not be detected.
-This is likely to discard many interesting genes, especially if the clusters are finely resolved with weak separation.
-To give a concrete example, consider a mixed population of CD4^+^-only, CD8^+^-only, double-positive and double-negative T cells.
-Neither _Cd4_ or _Cd8_ would be detected as subpopulation-specific markers because each gene is expressed in two subpopulations such that the minimum effect would be small.
-In practice, the minimum and maximum are most helpful for diagnosing discrepancies between the mean and median, rather than being used directly for ranking.
-
-Ranking genes by the min-rank is similiar in stringency to ranking by the maximum effect size, in that both will respond to strong DE in a single comparison.
-However, the min-rank is more useful as it ensures that a single comparison to another cluster with consistently large effects does not dominate the ranking.
-If we select all genes with min-ranks less than or equal to $T$, the resulting set is the union of the top $T$ genes from all pairwise comparisons.
-This guarantees that our set contains at least $T$ genes that can distinguish our cluster of interest from any other cluster,
-which permits a comprehensive determination of a cluster's identity.
-We demonstrate below for cluster `r target`, taking the top $T=5$ genes with the largest Cohen's $d$ from each comparison to display in Figure \@ref(fig:pbmc-grouped-ranked-heat).
-
-```{r pbmc-grouped-ranked-heat, fig.cap=sprintf("Heatmap of the centered average log-expression values for the top potential marker genes for cluster %s in the PBMC dataset. The set of markers was selected as those genes with Cohen's $d$-derived min-ranks less than or equal to 5.", target)}
-ordered <- chosen[order(chosen$rank.logFC.cohen),]
-top.ranked <- ordered[ordered$rank.logFC.cohen <= 5,]
-rownames(top.ranked)
-
-plotGroupedHeatmap(sce.pbmc, features=rownames(top.ranked), group="label",
- center=TRUE, zlim=c(-3, 3))
-```
-
-Our discussion above has focused mainly on potential markers that are upregulated in our cluster of interest, as these are the easiest to interpret and experimentally validate.
-However, it also means that any cluster defined by downregulation of a marker will not contain that gene among the top features.
-This is occasionally relevant for subtypes or other states that are defined by low expression of particular genes.
-In such cases, focusing on upregulation may yield a disappointing set of markers,
-and it may be worth examining some of the lowest-ranked genes to see if there is any consistent downregulation compared to other clusters.
-
-```{r}
-# Omitting the decreasing=TRUE to focus on negative effects.
-ordered <- chosen[order(chosen$median.logFC.cohen),1:4]
-head(ordered)
-```
-
-## Obtaining the full effects
-
-For more complex questions, we may need to interrogate effect sizes from specific comparisons of interest.
-To do so, we set `full.stats=TRUE` to obtain the effect sizes for all pairwise comparisons involving a particular cluster.
-This is returned in the form of a nested `DataFrame` for each effect size type -
-in the example below, `full.AUC` contains the AUCs for the comparisons between cluster `r target` and every other cluster.
-
-```{r}
-marker.info <- scoreMarkers(sce.pbmc, colLabels(sce.pbmc), full.stats=TRUE)
-chosen <- marker.info[["12"]]
-chosen$full.AUC
-```
-
-Say we want to identify the genes that distinguish cluster 12 from other clusters with high _LYZ_ expression.
-We subset `full.AUC` to the relevant comparisons and sort on our summary statistic of choice to obtain a ranking of markers within this subset.
-This allows us to easily characterize subtle differences between closely related clusters.
-To illustrate, we use the smallest rank from `computeMinRank()` to identify the top DE genes in cluster `r target` compared to the other _LYZ_-high clusters (Figure \@ref(fig:pbmc-grouped-ranked-heat2)).
-
-```{r pbmc-grouped-ranked-heat2, fig.cap=sprintf("Heatmap of the centered average log-expression values for the top potential marker genes for cluster %s relative to other _LYZ_-high clusters in the PBMC dataset. The set of markers was selected as those genes with AUC-derived min-ranks less than or equal to 10.", target)}
-lyz.high <- c("7", "8", "10", "12", "15", "18") # based on inspection of the previous Figure.
-subset <- chosen$full.AUC[,colnames(chosen$full.AUC) %in% lyz.high]
-to.show <- subset[computeMinRank(subset) <= 10,]
-to.show
-plotGroupedHeatmap(sce.pbmc[,colLabels(sce.pbmc) %in% lyz.high],
- features=rownames(to.show), group="label", center=TRUE, zlim=c(-3, 3))
-```
-
-```{r, echo=FALSE}
-# Checking that the others listed are, indeed, LYZ high.
-for (i in lyz.high) {
- out <- marker.info[[i]]
- stopifnot(out["LYZ", "median.AUC"] > 0.9)
-}
-```
-
-Similarly, we can use the full set of effect sizes to define our own summary statistic if the precomputed measures are too coarse.
-For example, we may be interested in markers that are upregulated against some percentage - say, 80% - of other clusters.
-This improves the cluster specificity of the ranking by being more stringent than the median yet not as stringent as the minimum.
-We achieve this by computing and sorting on the 20^th^ percentile of effect sizes, as shown below.
-
-```{r}
-stat <- rowQuantiles(as.matrix(chosen$full.AUC), p=0.2)
-chosen[order(stat, decreasing=TRUE), 1:4] # just showing the basic stats for brevity.
-```
-
-## Using a log-fold change threshold
-
-The Cohen's $d$ and AUC calculations consider both the magnitude of the difference between clusters as well as the variability within each cluster.
-If the variability is lower, it is possible for a gene to have a large effect size even if the magnitude of the difference is small.
-These genes tend to be somewhat uninformative for cell type identification despite their strong differential expression (e.g., ribosomal protein genes).
-We would prefer genes with larger log-fold changes between clusters, even if they have higher variability.
-
-To favor the detection of such genes, we can compute the effect sizes relative to a log-fold change threshold by setting `lfc=` in `scoreMarkers()`.
-The definition of Cohen's $d$ is generalized to the standardized difference between the observed log-fold change and the specified `lfc` threshold.
-Similarly, the AUC is redefined as the probability of randomly picking an expression value from one cluster that is greater than a random value from the other cluster plus `lfc`.
-A large positive Cohen's $d$ and an AUC above 0.5 can only be obtained if the observed log-fold change between clusters is significantly greater than `lfc`.
-We demonstrate below by obtaining the top markers for cluster 4 in the PBMC dataset with `lfc=2` (Figure \@ref(fig:pbmc-markers-4-lfc)).
-
-```{r pbmc-markers-4-lfc, fig.height=6, fig.width=6, fig.cap="Dot plot of the top potential marker genes (as determined by the mean AUC) for cluster 4 in the PBMC dataset. Each row corrresponds to a marker gene and each column corresponds to a cluster. The size of each dot represents the proportion of cells with detected expression of the gene in the cluster, while the color is proportional to the average expression across all cells in that cluster."}
-marker.info.lfc <- scoreMarkers(sce.pbmc, colLabels(sce.pbmc), lfc=2)
-chosen2 <- marker.info.lfc[["4"]] # another cluster for some variety.
-chosen2 <- chosen2[order(chosen2$mean.AUC, decreasing=TRUE),]
-chosen2[,c("self.average", "other.average", "mean.AUC")]
-plotDots(sce.pbmc, rownames(chosen2)[1:10], group="label")
-```
-
-Note that the interpretation of the AUC and Cohen's $d$ becomes slightly more complicated when `lfc` is non-zero.
-If `lfc` is positive, a positive Cohen's $d$ and an AUC above 0.5 represents upregulation.
-However, a negative Cohen's $d$ or AUC below 0.5 may not represent downregulation; it may just indicate that the observed log-fold change is less than the specified `lfc`.
-The converse applies when `lfc` is negative, where the only conclusive interpretation occurs for downregulated genes.
-For the most part, this complication is not too problematic for routine marker detection, as we are mostly interested in upregulated genes with large positive Cohen's $d$ and AUCs above 0.5.
-
-## Handling blocking factors {#marker-batch}
-
-Large studies may contain factors of variation that are known and not interesting (e.g., batch effects, sex differences).
-If these are not modelled, they can interfere with marker gene detection - most obviously by inflating the variance within each cluster, but also by distorting the log-fold changes if the cluster composition varies across levels of the blocking factor.
-To avoid these issues, we specify the blocking factor via the `block=` argument, as demonstrated below for the 416B data set.
-
-```{r, results='asis', echo=FALSE}
-extractFromPackage("lun-416b.Rmd", package="OSCA.workflows",
- chunk="clustering", objects="sce.416b")
-```
-
-```{r}
-m.out <- scoreMarkers(sce.416b, colLabels(sce.416b), block=sce.416b$block)
-```
-
-For each gene, each pairwise comparison between clusters is performed separately in each level of the blocking factor - in this case, the plate of origin.
-By comparing within each batch, we cancel out any batch effects so that they are not conflated with the biological differences between subpopulations.
-The effect sizes are then averaged across batches to obtain a single value per comparison, using a weighted mean that accounts for the number of cells involved in the comparison in each batch.
-A similar correction is applied to the mean log-expression and proportion of detected cells inside and outside each cluster.
-
-```{r blocked-markers-416b, fig.height=6, fig.width=6, fig.cap="Distribution of expression values across clusters for the top potential marker genes from cluster 1 in the 416B dataset. Each point represents a cell and is colored by the batch of origin."}
-demo <- m.out[["1"]]
-ordered <- demo[order(demo$median.logFC.cohen, decreasing=TRUE),]
-ordered[,1:4]
-
-plotExpression(sce.416b, features=rownames(ordered)[1:6],
- x="label", colour_by="block")
-```
-
-The `block=` argument works for all effect sizes shown above and is robust to differences in the log-fold changes or variance between batches.
-However, it assumes that each pair of clusters is present in at least one batch.
-In scenarios where cells from two clusters never co-occur in the same batch, the associated pairwise comparison will be impossible and is ignored during calculation of summary statistics.
-
-
-
-## Session Info {-}
-
-```{r sessionInfo, echo=FALSE, results='asis'}
-prettySessionInfo()
-```
diff --git a/inst/book/marker-detection.html b/inst/book/marker-detection.html
new file mode 100644
index 0000000..f24917e
--- /dev/null
+++ b/inst/book/marker-detection.html
@@ -0,0 +1,9 @@
+
+
+
+
+
+
+
+
+
diff --git a/inst/book/normalization.Rmd b/inst/book/normalization.Rmd
deleted file mode 100644
index 34eb9b6..0000000
--- a/inst/book/normalization.Rmd
+++ /dev/null
@@ -1,235 +0,0 @@
----
-output:
- html_document
-bibliography: ref.bib
----
-
-# Normalization
-
-```{r setup, echo=FALSE, results="asis"}
-library(rebook)
-chapterPreamble()
-```
-
-## Motivation
-
-Systematic differences in sequencing coverage between libraries are often observed in single-cell RNA sequencing data [@stegle2015computational].
-They typically arise from technical differences in cDNA capture or PCR amplification efficiency across cells, attributable to the difficulty of achieving consistent library preparation with minimal starting material.
-Normalization aims to remove these differences such that they do not interfere with comparisons of the expression profiles between cells.
-This ensures that any observed heterogeneity or differential expression within the cell population are driven by biology and not technical biases.
-
-We will mostly focus our attention on scaling normalization, which is the simplest and most commonly used class of normalization strategies.
-This involves dividing all counts for each cell by a cell-specific scaling factor, often called a "size factor" [@anders2010differential].
-The assumption here is that any cell-specific bias (e.g., in capture or amplification efficiency) affects all genes equally via scaling of the expected mean count for that cell.
-The size factor for each cell represents the estimate of the relative bias in that cell, so division of its counts by its size factor should remove that bias.
-The resulting "normalized expression values" can then be used for downstream analyses such as clustering and dimensionality reduction.
-To demonstrate, we will use the @zeisel2015brain dataset from the `r BiocStyle::Biocpkg("scRNAseq")` package.
-
-```{r, echo=FALSE, results="asis"}
-extractFromPackage("zeisel-brain.Rmd", package="OSCA.workflows",
- chunk="quality-control", objects="sce.zeisel")
-```
-
-```{r}
-sce.zeisel
-```
-
-## Library size normalization
-
-Library size normalization is the simplest strategy for performing scaling normalization.
-We define the library size as the total sum of counts across all genes for each cell, the expected value of which is assumed to scale with any cell-specific biases.
-The "library size factor" for each cell is then directly proportional to its library size where the proportionality constant is defined such that the mean size factor across all cells is equal to 1.
-This definition ensures that the normalized expression values are on the same scale as the original counts, which is useful for interpretation - especially when dealing with transformed data (see Section \@ref(normalization-transformation)).
-
-```{r}
-library(scater)
-lib.sf.zeisel <- librarySizeFactors(sce.zeisel)
-summary(lib.sf.zeisel)
-```
-
-In the Zeisel brain data, the library size factors differ by up to 10-fold across cells (Figure \@ref(fig:histlib)).
-This is typical of the variability in coverage in scRNA-seq data.
-
-```{r histlib, fig.cap="Distribution of size factors derived from the library size in the Zeisel brain dataset."}
-hist(log10(lib.sf.zeisel), xlab="Log10[Size factor]", col='grey80')
-```
-
-Strictly speaking, the use of library size factors assumes that there is no "imbalance" in the differentially expressed (DE) genes between any pair of cells.
-That is, any upregulation for a subset of genes is cancelled out by the same magnitude of downregulation in a different subset of genes.
-This ensures that the library size is an unbiased estimate of the relative cell-specific bias by avoiding composition effects [@robinson2010scaling].
-However, balanced DE is not generally present in scRNA-seq applications, which means that library size normalization may not yield accurate normalized expression values for downstream analyses.
-
-In practice, normalization accuracy is not a major consideration for exploratory scRNA-seq data analyses.
-Composition biases do not usually affect the separation of clusters, only the magnitude - and to a lesser extent, direction - of the log-fold changes between clusters or cell types.
-As such, library size normalization is usually sufficient in many applications where the aim is to identify clusters and the top markers that define each cluster.
-
-## Normalization by deconvolution
-
-As previously mentioned, composition biases will be present when any unbalanced differential expression exists between samples.
-Consider the simple example of two cells where a single gene $X$ is upregulated in one cell $A$ compared to the other cell $B$.
-This upregulation means that either (i) more sequencing resources are devoted to $X$ in $A$, thus decreasing coverage of all other non-DE genes when the total library size of each cell is experimentally fixed (e.g., due to library quantification);
-or (ii) the library size of $A$ increases when $X$ is assigned more reads or UMIs, increasing the library size factor and yielding smaller normalized expression values for all non-DE genes.
-In both cases, the net effect is that non-DE genes in $A$ will incorrectly appear to be downregulated compared to $B$.
-
-The removal of composition biases is a well-studied problem for bulk RNA sequencing data analysis.
-Normalization can be performed with the `estimateSizeFactorsFromMatrix()` function in the `r Biocpkg("DESeq2")` package [@anders2010differential;@love2014moderated] or with the `calcNormFactors()` function [@robinson2010scaling] in the `r Biocpkg("edgeR")` package.
-These assume that most genes are not DE between cells.
-Any systematic difference in count size across the non-DE majority of genes between two cells is assumed to represent bias that is used to compute an appropriate size factor for its removal.
-
-However, single-cell data can be problematic for these bulk normalization methods due to the dominance of low and zero counts.
-To overcome this, we pool counts from many cells to increase the size of the counts for accurate size factor estimation [@lun2016pooling].
-Pool-based size factors are then "deconvolved" into cell-based factors for normalization of each cell's expression profile.
-This is performed using the `calculateSumFactors()` function from `r Biocpkg("scran")`, as shown below.
-
-```{r}
-library(scran)
-set.seed(100)
-clust.zeisel <- quickCluster(sce.zeisel)
-table(clust.zeisel)
-
-deconv.sf.zeisel <- calculateSumFactors(sce.zeisel, cluster=clust.zeisel)
-summary(deconv.sf.zeisel)
-```
-
-We use a pre-clustering step with `quickCluster()` where cells in each cluster are normalized separately and the size factors are rescaled to be comparable across clusters.
-This avoids the assumption that most genes are non-DE across the entire population - only a non-DE majority is required between pairs of clusters, which is a weaker assumption for highly heterogeneous populations.
-By default, `quickCluster()` will use an approximate algorithm for PCA based on methods from the `r CRANpkg("irlba")` package.
-The approximation relies on stochastic initialization so we need to set the random seed (via `set.seed()`) for reproducibility.
-
-We see that the deconvolution size factors exhibit cell type-specific deviations from the library size factors in Figure \@ref(fig:deconv-zeisel).
-This is consistent with the presence of composition biases that are introduced by strong differential expression between cell types.
-Use of the deconvolution size factors adjusts for these biases to improve normalization accuracy for downstream applications.
-
-```{r deconv-zeisel, fig.cap="Deconvolution size factor for each cell in the Zeisel brain dataset, compared to the equivalent size factor derived from the library size. The red line corresponds to identity between the two size factors."}
-plot(lib.sf.zeisel, deconv.sf.zeisel, xlab="Library size factor",
- ylab="Deconvolution size factor", log='xy', pch=16,
- col=as.integer(factor(sce.zeisel$level1class)))
-abline(a=0, b=1, col="red")
-```
-
-Accurate normalization is most important for procedures that involve estimation and interpretation of per-gene statistics.
-For example, composition biases can compromise DE analyses by systematically shifting the log-fold changes in one direction or another.
-However, it tends to provide less benefit over simple library size normalization for cell-based analyses such as clustering.
-The presence of composition biases already implies strong differences in expression profiles, so changing the normalization strategy is unlikely to affect the outcome of a clustering procedure.
-
-## Normalization by spike-ins {#spike-norm}
-
-Spike-in normalization is based on the assumption that the same amount of spike-in RNA was added to each cell [@lun2017assessing].
-Systematic differences in the coverage of the spike-in transcripts can only be due to cell-specific biases, e.g., in capture efficiency or sequencing depth.
-To remove these biases, we equalize spike-in coverage across cells by scaling with "spike-in size factors".
-Compared to the previous methods, spike-in normalization requires no assumption about the biology of the system (i.e., the absence of many DE genes).
-Instead, it assumes that the spike-in transcripts were (i) added at a constant level to each cell, and (ii) respond to biases in the same relative manner as endogenous genes.
-
-Practically, spike-in normalization should be used if differences in the total RNA content of individual cells are of interest and must be preserved in downstream analyses.
-For a given cell, an increase in its overall amount of endogenous RNA will not increase its spike-in size factor.
-This ensures that the effects of total RNA content on expression across the population will not be removed upon scaling.
-By comparison, the other normalization methods described above will simply interpret any change in total RNA content as part of the bias and remove it.
-
-We demonstrate the use of spike-in normalization on a different dataset involving T cell activation after stimulation with T cell recepter ligands of varying affinity [@richard2018tcell].
-
-```{r}
-library(scRNAseq)
-sce.richard <- RichardTCellData()
-sce.richard <- sce.richard[,sce.richard$`single cell quality`=="OK"]
-sce.richard
-```
-
-We apply the `computeSpikeFactors()` method to estimate spike-in size factors for all cells.
-This is defined by converting the total spike-in count per cell into a size factor, using the same reasoning as in `librarySizeFactors()`.
-Scaling will subsequently remove any differences in spike-in coverage across cells.
-
-```{r}
-sce.richard <- computeSpikeFactors(sce.richard, "ERCC")
-summary(sizeFactors(sce.richard))
-```
-
-We observe a positive correlation between the spike-in size factors and deconvolution size factors within each treatment condition (Figure \@ref(fig:norm-spike-t)), indicating that they are capturing similar technical biases in sequencing depth and capture efficiency.
-However, we also observe that increasing stimulation of the T cell receptor - in terms of increasing affinity or time - results in a decrease in the spike-in factors relative to the library size factors.
-This is consistent with an increase in biosynthetic activity and total RNA content during stimulation, which reduces the relative spike-in coverage in each library (thereby decreasing the spike-in size factors) but increases the coverage of endogenous genes (thus increasing the library size factors).
-
-```{r norm-spike-t, fig.wide=TRUE, fig.asp=0.5, fig.width=8, fig.cap="Size factors from spike-in normalization, plotted against the library size factors for all cells in the T cell dataset. Each plot represents a different ligand treatment and each point is a cell coloured according by time from stimulation."}
-to.plot <- data.frame(
- DeconvFactor=calculateSumFactors(sce.richard),
- SpikeFactor=sizeFactors(sce.richard),
- Stimulus=sce.richard$stimulus,
- Time=sce.richard$time
-)
-
-ggplot(to.plot, aes(x=DeconvFactor, y=SpikeFactor, color=Time)) +
- geom_point() + facet_wrap(~Stimulus) + scale_x_log10() +
- scale_y_log10() + geom_abline(intercept=0, slope=1, color="red")
-```
-
-The differences between these two sets of size factors have real consequences for downstream interpretation.
-If the spike-in size factors were applied to the counts, the expression values in unstimulated cells would be scaled up while expression in stimulated cells would be scaled down.
-However, the opposite would occur if the deconvolution size factors were used.
-This can manifest as shifts in the magnitude and direction of DE between conditions when we switch between normalization strategies, as shown below for _Malat1_ (Figure \@ref(fig:norm-effect-malat)).
-
-```{r norm-effect-malat, fig.asp=1, fig.wide=TRUE, fig.cap="Distribution of log-normalized expression values for _Malat1_ after normalization with the deconvolution size factors (left) or spike-in size factors (right). Cells are stratified by the ligand affinity and colored by the time after stimulation."}
-# See below for explanation of logNormCounts().
-sce.richard.deconv <- logNormCounts(sce.richard, size_factors=to.plot$DeconvFactor)
-sce.richard.spike <- logNormCounts(sce.richard, size_factors=to.plot$SpikeFactor)
-
-gridExtra::grid.arrange(
- plotExpression(sce.richard.deconv, x="stimulus",
- colour_by="time", features="ENSMUSG00000092341") +
- theme(axis.text.x = element_text(angle = 90)) +
- ggtitle("After deconvolution"),
- plotExpression(sce.richard.spike, x="stimulus",
- colour_by="time", features="ENSMUSG00000092341") +
- theme(axis.text.x = element_text(angle = 90)) +
- ggtitle("After spike-in normalization"),
- ncol=2
-)
-```
-
-```{r, echo=FALSE, results="hide"}
-# Checking that Malat1 does indeed change in its DE status
-markers1 <- scoreMarkers(sce.richard.deconv, sce.richard$stimulus, full.stats=TRUE)
-markers2 <- scoreMarkers(sce.richard.spike, sce.richard$stimulus, full.stats=TRUE)
-stopifnot(
- all(as.matrix(markers1$unstimulated["ENSMUSG00000092341","full.logFC.cohen"]) -
- as.matrix(markers2$unstimulated["ENSMUSG00000092341","full.logFC.cohen"]) > 0.2)
-)
-```
-
-Whether or not total RNA content is relevant -- and thus, the choice of normalization strategy -- depends on the biological hypothesis.
-In most cases, changes in total RNA content are not interesting and can be normalized out by applying the library size or deconvolution factors.
-However, this may not always be appropriate if differences in total RNA are associated with a biological process of interest, e.g., cell cycle activity or T cell activation.
-Spike-in normalization will preserve these differences such that any changes in expression between biological groups have the correct sign.
-
-**However!**
-Regardless of whether we care about total RNA content, it is critical that the spike-in transcripts are normalized using the spike-in size factors.
-Size factors computed from the counts for endogenous genes should not be applied to the spike-in transcripts, precisely because the former captures differences in total RNA content that are not experienced by the latter.
-Attempting to normalize the spike-in counts with the gene-based size factors will lead to over-normalization and incorrect quantification.
-Thus, if normalized spike-in data is required, we must compute a separate set of size factors for the spike-in transcripts; this is automatically performed by functions such as `modelGeneVarWithSpikes()`.
-
-## Scaling and log-transforming {#normalization-transformation}
-
-Once we have computed the size factors, we use the `logNormCounts()` function from `r Biocpkg("scater")` to compute normalized expression values for each cell.
-This is done by dividing the count for each gene/spike-in transcript with the appropriate size factor for that cell.
-The function also log-transforms the normalized values, creating a new assay called `"logcounts"`.
-(Technically, these are "log-transformed normalized expression values", but that's too much of a mouthful to fit into the assay name.)
-These log-values will be the basis of our downstream analyses in the following chapters.
-
-```{r}
-set.seed(100)
-clust.zeisel <- quickCluster(sce.zeisel)
-sce.zeisel <- computeSumFactors(sce.zeisel, cluster=clust.zeisel, min.mean=0.1)
-sce.zeisel <- logNormCounts(sce.zeisel)
-assayNames(sce.zeisel)
-```
-
-The log-transformation is useful as differences in the log-values represent log-fold changes in expression.
-This is important in downstream procedures based on Euclidean distances, which includes many forms of clustering and dimensionality reduction.
-By operating on log-transformed data, we ensure that these procedures are measuring distances between cells based on log-fold changes in expression.
-Or in other words, which is more interesting - a gene that is expressed at an average count of 50 in cell type $A$ and 10 in cell type $B$, or a gene that is expressed at an average count of 1100 in $A$ and 1000 in $B$?
-Log-transformation focuses on the former by promoting contributions from genes with strong relative differences.
-
-See `r link("more-norm", "OSCA.advanced")` for further comments on transformation strategies.
-
-## Session Info {-}
-
-```{r sessionInfo, echo=FALSE, results='asis'}
-prettySessionInfo()
-```
diff --git a/inst/book/normalization.html b/inst/book/normalization.html
new file mode 100644
index 0000000..f3670be
--- /dev/null
+++ b/inst/book/normalization.html
@@ -0,0 +1,9 @@
+
+
+
+
+
+
+
+
+
diff --git a/inst/book/quality-control.Rmd b/inst/book/quality-control.Rmd
deleted file mode 100644
index e236406..0000000
--- a/inst/book/quality-control.Rmd
+++ /dev/null
@@ -1,292 +0,0 @@
----
-output:
- html_document
-bibliography: ref.bib
----
-
-# Quality Control
-
-```{r setup, echo=FALSE, results="asis"}
-library(rebook)
-chapterPreamble()
-```
-
-## Motivation {#quality-control-motivation}
-
-Low-quality libraries in scRNA-seq data can arise from a variety of sources such as cell damage during dissociation or failure in library preparation (e.g., inefficient reverse transcription or PCR amplification).
-These usually manifest as "cells" with low total counts, few expressed genes and high mitochondrial or spike-in proportions.
-These low-quality libraries are problematic as they can contribute to misleading results in downstream analyses:
-
-- They form their own distinct cluster(s), complicating interpretation of the results.
-This is most obviously driven by increased mitochondrial proportions or enrichment for nuclear RNAs after cell damage.
-In the worst case, low-quality libraries generated from different cell types can cluster together based on similarities in the damage-induced expression profiles, creating artificial intermediate states or trajectories between otherwise distinct subpopulations.
-Additionally, very small libraries can form their own clusters due to shifts in the mean upon transformation [@lun2018overcoming].
-- They interfere with characterization of population heterogeneity during variance estimation or principal components analysis.
-The first few principal components will capture differences in quality rather than biology, reducing the effectiveness of dimensionality reduction.
-Similarly, genes with the largest variances will be driven by differences between low- and high-quality cells.
-The most obvious example involves low-quality libraries with very low counts where scaling normalization inflates the apparent variance of genes that happen to have a non-zero count in those libraries.
-- They contain genes that appear to be strongly "upregulated" due to aggressive scaling to normalize for small library sizes.
-This is most problematic for contaminating transcripts (e.g., from the ambient solution) that are present in all libraries at low but constant levels.
-Increased scaling in low-quality libraries transforms small counts for these transcripts in large normalized expression values, resulting in apparent upregulation compared to other cells.
-This can be misleading as the affected genes are often biologically sensible but are actually expressed in another subpopulation.
-
-To avoid - or at least mitigate - these problems, we need to remove the problematic cells at the start of the analysis.
-This step is commonly referred to as quality control (QC) on the cells.
-(We will use "library" and "cell" rather interchangeably here, though the distinction will become important when dealing with droplet-based data.)
-We demonstrate using a small scRNA-seq dataset from @lun2017assessing, which is provided with no prior QC so that we can apply our own procedures.
-
-```{r, results='asis', echo=FALSE}
-extractFromPackage("lun-416b.Rmd", package="OSCA.workflows",
- chunk="loading", objects="sce.416b")
-```
-
-```{r}
-sce.416b
-```
-
-## Common choices of QC metrics
-
-We use several common QC metrics to identify low-quality cells based on their expression profiles.
-These metrics are described below in terms of reads for SMART-seq2 data, but the same definitions apply to UMI data generated by other technologies like MARS-seq and droplet-based protocols.
-
-- The library size is defined as the total sum of counts across all relevant features for each cell.
-Here, we will consider the relevant features to be the endogenous genes.
-Cells with small library sizes are of low quality as the RNA has been lost at some point during library preparation,
-either due to cell lysis or inefficient cDNA capture and amplification.
-- The number of expressed features in each cell is defined as the number of endogenous genes with non-zero counts for that cell.
-Any cell with very few expressed genes is likely to be of poor quality as the diverse transcript population has not been successfully captured.
-- The proportion of reads mapped to spike-in transcripts is calculated relative to the total count across all features (including spike-ins) for each cell.
-As the same amount of spike-in RNA should have been added to each cell, any enrichment in spike-in counts is symptomatic of loss of endogenous RNA.
-Thus, high proportions are indicative of poor-quality cells where endogenous RNA has been lost due to, e.g., partial cell lysis or RNA degradation during dissociation.
-- In the absence of spike-in transcripts, the proportion of reads mapped to genes in the mitochondrial genome can be used.
-High proportions are indicative of poor-quality cells [@islam2014quantitative;@ilicic2016classification], presumably because of loss of cytoplasmic RNA from perforated cells.
-The reasoning is that, in the presence of modest damage, the holes in the cell membrane permit efflux of individual transcript molecules but are too small to allow mitochondria to escape, leading to a relative enrichment of mitochondrial transcripts.
-For single-nuclei RNA-seq experiments, high proportions are also useful as they can mark cells where the cytoplasm has not been successfully stripped.
-
-For each cell, we calculate these QC metrics using the `perCellQCMetrics()` function from the `r BiocStyle::Biocpkg("scater")` package [@mccarthy2017scater].
-The `sum` column contains the total count for each cell and the `detected` column contains the number of detected genes.
-The `subsets_Mito_percent` column contains the percentage of reads mapped to mitochondrial transcripts.
-Finally, the `altexps_ERCC_percent` column contains the percentage of reads mapped to ERCC transcripts.
-
-```{r}
-library(SingleCellExperiment)
-# Identifying the mitochondrial transcripts in our SingleCellExperiment.
-location <- rowRanges(sce.416b)
-is.mito <- any(seqnames(location)=="MT")
-
-library(scuttle)
-df <- perCellQCMetrics(sce.416b, subsets=list(Mito=is.mito))
-summary(df$sum)
-summary(df$detected)
-summary(df$subsets_Mito_percent)
-summary(df$altexps_ERCC_percent)
-```
-
-Alternatively, users may prefer to use the `addPerCellQC()` function.
-This computes and appends the per-cell QC statistics to the `colData` of the `SingleCellExperiment` object,
-allowing us to retain all relevant information in a single object for later manipulation.
-
-```{r}
-sce.416b <- addPerCellQCMetrics(sce.416b, subsets=list(Mito=is.mito))
-colnames(colData(sce.416b))
-```
-
-A key assumption here is that the QC metrics are independent of the biological state of each cell.
-Poor values (e.g., low library sizes, high mitochondrial proportions) are presumed to be driven by technical factors rather than biological processes, meaning that the subsequent removal of cells will not misrepresent the biology in downstream analyses.
-Major violations of this assumption would potentially result in the loss of cell types that have, say, systematically low RNA content or high numbers of mitochondria.
-We can check for such violations using diagnostic plots described in Section \@ref(quality-control-plots) and `r link("qc-discard-cell-types", "OSCA.advanced")`.
-
-## Identifying low-quality cells
-
-### With fixed thresholds {#fixed-qc}
-
-The simplest approach to identifying low-quality cells involves applying fixed thresholds to the QC metrics.
-For example, we might consider cells to be low quality if they have library sizes below 100,000 reads; express fewer than 5,000 genes; have spike-in proportions above 10%; or have mitochondrial proportions above 10%.
-
-```{r}
-qc.lib <- df$sum < 1e5
-qc.nexprs <- df$detected < 5e3
-qc.spike <- df$altexps_ERCC_percent > 10
-qc.mito <- df$subsets_Mito_percent > 10
-discard <- qc.lib | qc.nexprs | qc.spike | qc.mito
-
-# Summarize the number of cells removed for each reason.
-DataFrame(LibSize=sum(qc.lib), NExprs=sum(qc.nexprs),
- SpikeProp=sum(qc.spike), MitoProp=sum(qc.mito), Total=sum(discard))
-```
-
-While simple, this strategy requires considerable experience to determine appropriate thresholds for each experimental protocol and biological system.
-Thresholds for read count-based data are not applicable for UMI-based data, and vice versa.
-Differences in mitochondrial activity or total RNA content require constant adjustment of the mitochondrial and spike-in thresholds, respectively, for different biological systems.
-Indeed, even with the same protocol and system, the appropriate threshold can vary from run to run due to the vagaries of cDNA capture efficiency and sequencing depth per cell.
-
-### With adaptive thresholds {#quality-control-outlier}
-
-Here, we assume that most of the dataset consists of high-quality cells.
-We then identify cells that are outliers for the various QC metrics, based on the median absolute deviation (MAD) from the median value of each metric across all cells.
-By default, we consider a value to be an outlier if it is more than 3 MADs from the median in the "problematic" direction.
-This is loosely motivated by the fact that such a filter will retain 99% of non-outlier values that follow a normal distribution.
-We demonstrate by using the `perCellQCFilters()` function on the QC metrics from the 416B dataset.
-
-```{r}
-reasons <- perCellQCFilters(df,
- sub.fields=c("subsets_Mito_percent", "altexps_ERCC_percent"))
-colSums(as.matrix(reasons))
-```
-
-This function will identify cells with log-transformed library sizes that are more than 3 MADs below the median.
-A log-transformation is used to improve resolution at small values when `type="lower"`
-and to avoid negative thresholds that would be meaningless for a non-negative metric.
-Furthermore, it is not uncommon for the distribution of library sizes to exhibit a heavy right tail;
-the log-transformation avoids inflation of the MAD in a manner that might compromise outlier detection on the left tail.
-(More generally, it makes the distribution seem more normal to justify the 99% rationale mentioned above.)
-The function will also do the same for the log-transformed number of expressed genes.
-
-`perCellQCFilters()` will also identify outliers for the proportion-based metrics specified in the `sub.fields=` arguments.
-These distributions frequently exhibit a heavy right tail, but unlike the two previous metrics, it is the right tail itself that contains the putative low-quality cells.
-Thus, we do not perform any transformation to shrink the tail - rather, our hope is that the cells in the tail are identified as large outliers.
-(While it is theoretically possible to obtain a meaningless threshold above 100%, this is rare enough to not be of practical concern.)
-
-A cell that is an outlier for any of these metrics is considered to be of low quality and discarded.
-This is captured in the `discard` column, which can be used for later filtering (Section \@ref(quality-control-discarded)).
-
-```{r}
-summary(reasons$discard)
-```
-
-We can also extract the exact filter thresholds from the attributes of each of the logical vectors.
-This may be useful for checking whether the automatically selected thresholds are appropriate.
-
-```{r}
-attr(reasons$low_lib_size, "thresholds")
-attr(reasons$low_n_features, "thresholds")
-```
-
-With this strategy, the thresholds adapt to both the location and spread of the distribution of values for a given metric.
-This allows the QC procedure to adjust to changes in sequencing depth, cDNA capture efficiency, mitochondrial content, etc. without requiring any user intervention or prior experience.
-However, the underlying assumption of a high-quality majority may not always be appropriate, which is discussed in more detail in `r link("outlier-assumptions", "OSCA.advanced")`.
-
-### Other approaches
-
-Another strategy is to identify outliers in high-dimensional space based on the QC metrics for each cell.
-We use methods from `r CRANpkg("robustbase")` to quantify the "outlyingness" of each cells based on their QC metrics, and then use `isOutlier()` to identify low-quality cells that exhibit unusually high levels of outlyingness.
-
-```{r}
-stats <- cbind(log10(df$sum), log10(df$detected),
- df$subsets_Mito_percent, df$altexps_ERCC_percent)
-
-library(robustbase)
-outlying <- adjOutlyingness(stats, only.outlyingness = TRUE)
-multi.outlier <- isOutlier(outlying, type = "higher")
-summary(multi.outlier)
-```
-
-This and related approaches like PCA-based outlier detection and support vector machines can provide more power to distinguish low-quality cells from high-quality counterparts [@ilicic2016classification] as they can exploit patterns across many QC metrics.
-However, this comes at some cost to interpretability, as the reason for removing a given cell may not always be obvious.
-
-For completeness, we note that outliers can also be identified from the gene expression profiles, rather than QC metrics.
-We consider this to be a risky strategy as it can remove high-quality cells in rare populations.
-
-## Checking diagnostic plots {#quality-control-plots}
-
-It is good practice to inspect the distributions of QC metrics (Figure \@ref(fig:qc-dist-416b)) to identify possible problems.
-In the most ideal case, we would see normal distributions that would justify the 3 MAD threshold used in outlier detection.
-A large proportion of cells in another mode suggests that the QC metrics might be correlated with some biological state, potentially leading to the loss of distinct cell types during filtering;
-or that there were inconsistencies with library preparation for a subset of cells, a not-uncommon phenomenon in plate-based protocols.
-
-```{r qc-dist-416b, fig.asp=2, fig.wide=TRUE, fig.cap="Distribution of QC metrics for each batch and phenotype in the 416B dataset. Each point represents a cell and is colored according to whether it was discarded, respectively."}
-colData(sce.416b) <- cbind(colData(sce.416b), df)
-sce.416b$block <- factor(sce.416b$block)
-sce.416b$phenotype <- ifelse(grepl("induced", sce.416b$phenotype),
- "induced", "wild type")
-sce.416b$discard <- reasons$discard
-
-library(scater)
-gridExtra::grid.arrange(
- plotColData(sce.416b, x="block", y="sum", colour_by="discard",
- other_fields="phenotype") + facet_wrap(~phenotype) +
- scale_y_log10() + ggtitle("Total count"),
- plotColData(sce.416b, x="block", y="detected", colour_by="discard",
- other_fields="phenotype") + facet_wrap(~phenotype) +
- scale_y_log10() + ggtitle("Detected features"),
- plotColData(sce.416b, x="block", y="subsets_Mito_percent",
- colour_by="discard", other_fields="phenotype") +
- facet_wrap(~phenotype) + ggtitle("Mito percent"),
- plotColData(sce.416b, x="block", y="altexps_ERCC_percent",
- colour_by="discard", other_fields="phenotype") +
- facet_wrap(~phenotype) + ggtitle("ERCC percent"),
- ncol=1
-)
-```
-
-Another useful diagnostic involves plotting the proportion of mitochondrial counts against some of the other QC metrics.
-The aim is to confirm that there are no cells with both large total counts and large mitochondrial counts, to ensure that we are not inadvertently removing high-quality cells that happen to be highly metabolically active (e.g., hepatocytes).
-We demonstrate using data from a larger experiment involving the mouse brain [@zeisel2015brain];
-in this case, we do not observe any points in the top-right corner in Figure \@ref(fig:qc-mito-zeisel) that might potentially correspond to metabolically active, undamaged cells.
-
-```{r, echo=FALSE, results="asis"}
-extractFromPackage("zeisel-brain.Rmd", package="OSCA.workflows",
- chunk="gene-annotation", objects="sce.zeisel")
-```
-
-```{r qc-mito-zeisel, fig.cap="Percentage of UMIs assigned to mitochondrial transcripts in the Zeisel brain dataset, plotted against the total number of UMIs (top). Each point represents a cell and is colored according to whether it was considered low-quality and discarded."}
-sce.zeisel <- addPerCellQC(sce.zeisel,
- subsets=list(Mt=rowData(sce.zeisel)$featureType=="mito"))
-
-qc <- quickPerCellQC(colData(sce.zeisel),
- sub.fields=c("altexps_ERCC_percent", "subsets_Mt_percent"))
-sce.zeisel$discard <- qc$discard
-
-plotColData(sce.zeisel, x="sum", y="subsets_Mt_percent", colour_by="discard")
-```
-
-Comparison of the ERCC and mitochondrial percentages can also be informative (Figure \@ref(fig:qc-mito-spike-zeisel)).
-Low-quality cells with small mitochondrial percentages, large spike-in percentages and small library sizes are likely to be stripped nuclei, i.e., they have been so extensively damaged that they have lost all cytoplasmic content.
-On the other hand, cells with high mitochondrial percentages and low ERCC percentages may represent undamaged cells that are metabolically active.
-This interpretation also applies for single-nuclei studies but with a switch of focus:
-the stripped nuclei become the libraries of interest while the undamaged cells are considered to be low quality.
-
-```{r qc-mito-spike-zeisel, fig.cap="Percentage of UMIs assigned to mitochondrial transcripts in the Zeisel brain dataset, plotted against the percentage of UMIs assigned to spike-in transcripts (bottom). Each point represents a cell and is colored according to whether it was considered low-quality and discarded."}
-plotColData(sce.zeisel, x="altexps_ERCC_percent", y="subsets_Mt_percent",
- colour_by="discard")
-```
-
-We see that all of these metrics exhibit weak correlations to each other,
-presumably a manifestation of a common underlying effect of cell damage.
-The weakness of the correlations motivates the use of several metrics to capture different aspects of technical quality.
-Of course, the flipside is that these metrics may also represent different aspects of biology, increasing the risk of inadvertently discarding entire cell types.
-
-## Removing low-quality cells {#quality-control-discarded}
-
-Once low-quality cells have been identified, we can choose to either remove them or mark them.
-Removal is the most straightforward option and is achieved by subsetting the `SingleCellExperiment` by column.
-In this case, we use the low-quality calls from Section \@ref(quality-control-outlier) to generate a subsetted `SingleCellExperiment` that we would use for downstream analyses.
-
-```{r}
-# Keeping the columns we DON'T want to discard.
-filtered <- sce.416b[,!reasons$discard]
-```
-
-The other option is to simply mark the low-quality cells as such and retain them in the downstream analysis.
-The aim here is to allow clusters of low-quality cells to form, and then to identify and ignore such clusters during interpretation of the results.
-This approach avoids discarding cell types that have poor values for the QC metrics, deferring the decision on whether a cluster of such cells represents a genuine biological state.
-
-```{r}
-marked <- sce.416b
-marked$discard <- reasons$discard
-```
-
-The downside is that it shifts the burden of QC to the manual interpretation of the clusters, which is already a major bottleneck in scRNA-seq data analysis (Chapters \@ref(clustering), \@ref(marker-detection) and \@ref(cell-type-annotation)).
-Indeed, if we do not trust the QC metrics, we would have to distinguish between genuine cell types and low-quality cells based only on marker genes, and this is not always easy due to the tendency of the latter to "express" interesting genes (Section \@ref(quality-control-motivation)).
-Retention of low-quality cells also compromises the accuracy of the variance modelling, requiring, e.g., use of more PCs to offset the fact that the early PCs are driven by differences between low-quality and other cells.
-
-For routine analyses, we suggest performing removal by default to avoid complications from low-quality cells.
-This allows most of the population structure to be characterized with no - or, at least, fewer - concerns about its validity.
-Once the initial analysis is done, and if there are any concerns about discarded cell types (`r link("qc-discard-cell-types", "OSCA.advanced")`), a more thorough re-analysis can be performed where the low-quality cells are only marked.
-This recovers cell types with low RNA content, high mitochondrial proportions, etc. that only need to be interpreted insofar as they "fill the gaps" in the initial analysis.
-
-## Session Info {-}
-
-```{r sessionInfo, echo=FALSE, results='asis'}
-prettySessionInfo()
-```
diff --git a/inst/book/quality-control.html b/inst/book/quality-control.html
new file mode 100644
index 0000000..c842b1e
--- /dev/null
+++ b/inst/book/quality-control.html
@@ -0,0 +1,9 @@
+
+
+
+
+
+
+
+
+
diff --git a/inst/book/reduced-dimensions.Rmd b/inst/book/reduced-dimensions.Rmd
deleted file mode 100644
index 6e7ddee..0000000
--- a/inst/book/reduced-dimensions.Rmd
+++ /dev/null
@@ -1,267 +0,0 @@
----
-output:
- html_document
-bibliography: ref.bib
----
-
-# Dimensionality reduction
-
-```{r setup, echo=FALSE, results="asis"}
-library(rebook)
-chapterPreamble()
-```
-
-## Overview
-
-Many scRNA-seq analysis procedures involve comparing cells based on their expression values across multiple genes.
-For example, clustering aims to identify cells with similar transcriptomic profiles by computing Euclidean distances across genes.
-In these applications, each individual gene represents a dimension of the data.
-More intuitively, if we had a scRNA-seq data set with two genes, we could make a two-dimensional plot where each axis represents the expression of one gene and each point in the plot represents a cell.
-This concept can be extended to data sets with thousands of genes where each cell's expression profile defines its location in the high-dimensional expression space.
-
-As the name suggests, dimensionality reduction aims to reduce the number of separate dimensions in the data.
-This is possible because different genes are correlated if they are affected by the same biological process.
-Thus, we do not need to store separate information for individual genes, but can instead compress multiple features into a single dimension, e.g., an "eigengene" [@langfelder2007eigengene].
-This reduces computational work in downstream analyses like clustering, as calculations only need to be performed for a few dimensions rather than thousands of genes;
-reduces noise by averaging across multiple genes to obtain a more precise representation of the patterns in the data;
-and enables effective plotting of the data, for those of us who are not capable of visualizing more than 3 dimensions.
-
-We will use the @zeisel2015brain dataset to demonstrate the applications of various dimensionality reduction methods in this chapter.
-
-```{r, echo=FALSE, results="asis"}
-extractFromPackage("zeisel-brain.Rmd", package="OSCA.workflows",
- chunk="variance-modelling",
- objects=c("sce.zeisel", "dec.zeisel", "top.hvgs"))
-```
-
-```{r}
-sce.zeisel
-```
-
-## Principal components analysis
-
-Principal components analysis (PCA) discovers axes in high-dimensional space that capture the largest amount of variation.
-This is best understood by imagining each axis as a line.
-Say we draw a line anywhere, and we move each cell in our data set onto the closest position on the line.
-The variance captured by this axis is defined as the variance in the positions of cells along that line.
-In PCA, the first axis (or "principal component", PC) is chosen such that it maximizes this variance.
-The next PC is chosen such that it is orthogonal to the first and captures the greatest remaining amount of variation, and so on.
-
-By definition, the top PCs capture the dominant factors of heterogeneity in the data set.
-In the context of scRNA-seq, our assumption is that biological processes affect multiple genes in a coordinated manner.
-This means that the earlier PCs are likely to represent biological structure as more variation can be captured by considering the correlated behavior of many genes.
-By comparison, random technical or biological noise is expected to affect each gene independently.
-There is unlikely to be an axis that can capture random variation across many genes, meaning that noise should mostly be concentrated in the later PCs.
-This motivates the use of the earlier PCs in our downstream analyses, which concentrates the biological signal to simultaneously reduce computational work and remove noise.
-
-The use of the earlier PCs for denoising and data compaction is a strategy that is simple, highly effective and widely used in a variety of fields.
-It takes advantage of the well-studied theoretical properties of the PCA -
-namely, that a low-rank approximation formed from the top PCs is the optimal approximation of the original data for a given matrix rank.
-Indeed, the Euclidean distances between cells in PC space can be treated as an approximation of the same distances in the original dataset
-The literature for PCA also provides us with a range of fast implementations for scalable and efficient data analysis.
-
-We perform a PCA on the log-normalized expression values using the `fixedPCA()` function from `r Biocpkg("scran")`.
-By default, `fixedPCA()` will compute the first 50 PCs and store them in the `reducedDims()` of the output `SingleCellExperiment` object, as shown below.
-Here, we use only the top 2000 genes with the largest biological components to reduce both computational work and high-dimensional random noise.
-In particular, while PCA is robust to random noise, an excess of it may cause the earlier PCs to capture noise instead of biological structure [@johnstone2009consistency].
-This effect can be mitigated by restricting the PCA to a subset of HVGs, for which we can use any of the strategies described in Chapter \@ref(feature-selection).
-
-```{r}
-library(scran)
-top.zeisel <- getTopHVGs(dec.zeisel, n=2000)
-
-set.seed(100) # See below.
-sce.zeisel <- fixedPCA(sce.zeisel, subset.row=top.zeisel)
-reducedDimNames(sce.zeisel)
-dim(reducedDim(sce.zeisel, "PCA"))
-```
-
-For large data sets, greater efficiency is obtained by using approximate SVD algorithms that only compute the top PCs.
-By default, most PCA-related functions in `r Biocpkg("scater")` and `r Biocpkg("scran")` will use methods from the `r CRANpkg("irlba")` or `r CRANpkg("rsvd")` packages to perform the SVD.
-We can explicitly specify the SVD algorithm to use by passing an `BiocSingularParam` object to the `BSPARAM=` argument (see `r link("big-data-svd", "OSCA.advanced")` for more details).
-Many of these approximate algorithms are based on randomization and thus require `set.seed()` to obtain reproducible results.
-
-```{r}
-library(BiocSingular)
-set.seed(1000)
-sce.zeisel <- fixedPCA(sce.zeisel, subset.row=top.zeisel,
- BSPARAM=RandomParam(), name="randomized")
-reducedDimNames(sce.zeisel)
-```
-
-## Choosing the number of PCs
-
-How many of the top PCs should we retain for downstream analyses?
-The choice of the number of PCs $d$ is a decision that is analogous to the choice of the number of HVGs to use.
-Using more PCs will retain more biological signal at the cost of including more noise that might mask said signal.
-On the other hand, using fewer PCs will introduce competition between different factors of variation,
-where weaker (but still interesting) factors may be pushed down into lower PCs and inadvertently discarded from downtream analyses.
-
-Much like the choice of the number of HVGs, it is hard to determine whether an "optimal" choice exists for the number of PCs.
-Certainly, we could attempt to remove the technical variation that is almost always uninteresting.
-However, even if we were only left with biological variation, there is no straightforward way to automatically determine which aspects of this variation are relevant.
-One analyst's biological signal may be irrelevant noise to another analyst with a different scientific question.
-For example, heterogeneity within a population might be interesting when studying continuous processes like metabolic flux or differentiation potential,
-but is comparable to noise in applications that only aim to distinguish between distinct cell types.
-
-Most practitioners will simply set $d$ to a "reasonable" but arbitrary value, typically ranging from 10 to 50.
-This is often satisfactory as the later PCs explain so little variance that their inclusion or omission has no major effect.
-For example, in the Zeisel dataset, few PCs explain more than 1\% of the variance in the entire dataset (Figure \@ref(fig:zeisel-scree)) and choosing between, say, 20 and 40 PCs would not even amount to four percentage points' worth of difference in variance.
-In fact, the main consequence of using more PCs is simply that downstream calculations take longer as they need to compute over more dimensions, but most PC-related calculations are fast enough that this is not a practical concern.
-
-```{r zeisel-scree, fig.cap="Percentage of variance explained by successive PCs in the Zeisel dataset, shown on a log-scale for visualization purposes."}
-percent.var <- attr(reducedDim(sce.zeisel), "percentVar")
-plot(percent.var, log="y", xlab="PC", ylab="Variance explained (%)")
-```
-
-```{r, echo=FALSE}
-stopifnot(sum(percent.var[20:40]) < 4)
-```
-
-Nonetheless, `r link("more-choices-for-the-number-of-pcs", "OSCA.advanced")` describes some more data-driven strategies to guide a suitable choice of $d$.
-These automated choices are best treated as guidelines as they make assumptions about what variation is "interesting".
-Indeed, the concepts in `r link("based-on-population-structure", "OSCA.advanced")` could even be used to provide some justification for an arbitrarily chosen $d$.
-More diligent readers may consider repeating the analysis with a variety of choices of $d$ to explore other perspectives of the dataset at a different bias-variance trade-off,
-though this tends to be unnecessary work in most applications.
-
-## Visualizing the PCs
-
-Algorithms are more than happy to operate on 10-50 PCs, but these are still too many dimensions for human comprehension.
-To visualize the data, we could take the top 2 PCs for plotting (Figure \@ref(fig:zeisel-pca)).
-
-```{r zeisel-pca, fig.cap="PCA plot of the first two PCs in the Zeisel brain data. Each point is a cell, coloured according to the annotation provided by the original authors."}
-library(scater)
-plotReducedDim(sce.zeisel, dimred="PCA", colour_by="level1class")
-```
-
-The problem is that PCA is a linear technique, i.e., only variation along a line in high-dimensional space is captured by each PC.
-As such, it cannot efficiently pack differences in $d$ dimensions into the first 2 PCs.
-This is demonstrated in Figure \@ref(fig:zeisel-pca) where the top two PCs fail to resolve some subpopulations identified by @zeisel2015brain.
-If the first PC is devoted to resolving the biggest difference between subpopulations, and the second PC is devoted to resolving the next biggest difference, then the remaining differences will not be visible in the plot.
-
-One workaround is to plot several of the top PCs against each other in pairwise plots (Figure \@ref(fig:zeisel-pca-multi)).
-However, it is difficult to interpret multiple plots simultaneously, and even this approach is not sufficient to separate some of the annotated subpopulations.
-
-```{r zeisel-pca-multi, fig.wide=TRUE, fig.asp=1, fig.cap="PCA plot of the first two PCs in the Zeisel brain data. Each point is a cell, coloured according to the annotation provided by the original authors."}
-plotReducedDim(sce.zeisel, dimred="PCA", ncomponents=4,
- colour_by="level1class")
-```
-
-Thus, plotting the top few PCs is not satisfactory for visualization of complex populations.
-That said, the PCA itself is still of great value in visualization as it compacts and denoises the data prior to downstream steps.
-The top PCs are often used as input to more sophisticated (and computationally intensive) algorithms for dimensionality reduction.
-
-## Non-linear methods for visualization
-
-### $t$-stochastic neighbor embedding
-
-The _de facto_ standard for visualization of scRNA-seq data is the $t$-stochastic neighbor embedding ($t$-SNE) method [@van2008visualizing].
-This attempts to find a low-dimensional representation of the data that preserves the distances between each point and its neighbors in the high-dimensional space.
-Unlike PCA, it is not restricted to linear transformations, nor is it obliged to accurately represent distances between distant populations.
-This means that it has much more freedom in how it arranges cells in low-dimensional space, enabling it to separate many distinct clusters in a complex population (Figure \@ref(fig:tsne-brain)).
-
-```{r tsne-brain, fig.cap="$t$-SNE plots constructed from the top PCs in the Zeisel brain dataset. Each point represents a cell, coloured according to the published annotation."}
-set.seed(00101001101)
-
-# runTSNE() stores the t-SNE coordinates in the reducedDims
-# for re-use across multiple plotReducedDim() calls.
-sce.zeisel <- runTSNE(sce.zeisel, dimred="PCA")
-plotReducedDim(sce.zeisel, dimred="TSNE", colour_by="level1class")
-```
-
-One of the main disadvantages of $t$-SNE is that it is much more computationally intensive than other visualization methods.
-We mitigate this effect by setting `dimred="PCA"` in `runTSNE()`, which instructs the function to perform the $t$-SNE calculations on the top PCs in `sce.zeisel`.
-This exploits the data compaction and noise removal of the PCA for faster and cleaner results in the $t$-SNE.
-It is also possible to run $t$-SNE on the original expression matrix but this is less efficient.
-
-Another issue with $t$-SNE is that it requires the user to be aware of additional parameters (discussed [here](http://distill.pub/2016/misread-tsne/) in some depth).
-It involves a random initialization so we need to set the seed to ensure that the chosen results are reproducible.
-We may also wish to repeat the visualization several times to ensure that the results are representative.
-
-The "perplexity" is another important parameter that determines the granularity of the visualization (Figure \@ref(fig:tsne-perplexity)).
-Low perplexities will favor resolution of finer structure, possibly to the point that the visualization is compromised by random noise.
-Thus, it is advisable to test different perplexity values to ensure that the choice of perplexity does not drive the interpretation of the plot.
-
-```{r tsne-perplexity, fig.cap="$t$-SNE plots constructed from the top PCs in the Zeisel brain dataset, using a range of perplexity values. Each point represents a cell, coloured according to its annotation.", fig.width=15, fig.asp=0.3}
-set.seed(100)
-sce.zeisel <- runTSNE(sce.zeisel, dimred="PCA", perplexity=5)
-out5 <- plotReducedDim(sce.zeisel, dimred="TSNE",
- colour_by="level1class") + ggtitle("perplexity = 5")
-
-set.seed(100)
-sce.zeisel <- runTSNE(sce.zeisel, dimred="PCA", perplexity=20)
-out20 <- plotReducedDim(sce.zeisel, dimred="TSNE",
- colour_by="level1class") + ggtitle("perplexity = 20")
-
-set.seed(100)
-sce.zeisel <- runTSNE(sce.zeisel, dimred="PCA", perplexity=80)
-out80 <- plotReducedDim(sce.zeisel, dimred="TSNE",
- colour_by="level1class") + ggtitle("perplexity = 80")
-
-gridExtra::grid.arrange(out5, out20, out80, ncol=3)
-```
-
-Finally, it is unwise to read too much into the relative sizes and positions of the visual clusters.
-$t$-SNE will inflate dense clusters and compress sparse ones, such that we cannot use the size as a measure of subpopulation heterogeneity.
-In addition, $t$-SNE is not obliged to preserve the relative locations of non-neighboring clusters, such that we cannot use their positions to determine relationships between distant clusters.
-
-Despite its shortcomings, $t$-SNE is proven tool for general-purpose visualization of scRNA-seq data and remains a popular choice in many analysis pipelines.
-In particular, this author enjoys looking at $t$-SNEs as they remind him of histology slides, which allows him to pretend that he is looking at real data.
-
-### Uniform manifold approximation and projection
-
-The uniform manifold approximation and projection (UMAP) method [@mcInnes2018umap] is an alternative to $t$-SNE for non-linear dimensionality reduction.
-It is roughly similar to $t$-SNE in that it also tries to find a low-dimensional representation that preserves relationships between neighbors in high-dimensional space.
-However, the two methods are based on different theory, represented by differences in the various graph weighting equations.
-This manifests as a different visualization as shown in Figure \@ref(fig:umap-brain).
-
-```{r umap-brain, fig.cap="UMAP plots constructed from the top PCs in the Zeisel brain dataset. Each point represents a cell, coloured according to the published annotation."}
-set.seed(1100101001)
-sce.zeisel <- runUMAP(sce.zeisel, dimred="PCA")
-plotReducedDim(sce.zeisel, dimred="UMAP", colour_by="level1class")
-```
-
-Compared to $t$-SNE, the UMAP visualization tends to have more compact visual clusters with more empty space between them.
-It also attempts to preserve more of the global structure than $t$-SNE.
-From a practical perspective, UMAP is much faster than $t$-SNE, which may be an important consideration for large datasets.
-(Nonetheless, we have still run UMAP on the top PCs here for consistency.)
-UMAP also involves a series of randomization steps so setting the seed is critical.
-
-Like $t$-SNE, UMAP has its own suite of hyperparameters that affect the visualization (see the documentation [here](https://umap-learn.readthedocs.io/en/latest/parameters.html)).
-Of these, the number of neighbors (`n_neighbors`) and the minimum distance between embedded points (`min_dist`) have the greatest effect on the granularity of the output.
-If these values are too low, random noise will be incorrectly treated as high-resolution structure, while values that are too high will discard fine structure altogether in favor of obtaining an accurate overview of the entire dataset.
-Again, it is a good idea to test a range of values for these parameters to ensure that they do not compromise any conclusions drawn from a UMAP plot.
-
-It is arguable whether the UMAP or $t$-SNE visualizations are more useful or aesthetically pleasing.
-UMAP aims to preserve more global structure but this necessarily reduces resolution within each visual cluster.
-However, UMAP is unarguably much faster, and for that reason alone, it is increasingly displacing $t$-SNE as the method of choice for visualizing large scRNA-seq data sets.
-
-### Interpreting the plots {#visualization-interpretation}
-
-Dimensionality reduction for visualization necessarily involves discarding information and distorting the distances between cells to fit high-dimensional data into a 2-dimensional space.
-One might wonder whether the results of such extreme data compression can be trusted.
-Indeed, some of our more quantitative colleagues consider such visualizations to be more artistic than scientific, fit for little but impressing collaborators and reviewers!
-Perhaps this perspective is not entirely invalid, but we suggest that there is some value to be extracted from them provided that they are accompanied by an analysis of a higher-rank representation.
-
-As a general rule, focusing on local neighborhoods provides the safest interpretation of $t$-SNE and UMAP plots.
-These methods spend considerable effort to ensure that each cell's nearest neighbors in the input high-dimensional space are still its neighbors in the output two-dimensional embedding.
-Thus, if we see multiple cell types or clusters in a single unbroken "island" in the embedding, we could infer that those populations were also close neighbors in higher-dimensional space.
-However, less can be said about the distances between non-neighboring cells; there is no guarantee that large distances are faithfully recapitulated in the embedding, given the distortions necessary for this type of dimensionality reduction.
-It would be courageous to use the distances between islands (seen to be measured, on occasion, with a ruler!) to make statements about the relative similarity of distinct cell types.
-
-On a related note, we prefer to restrict the $t$-SNE/UMAP coordinates for visualization and use the higher-rank representation for any quantitative analyses.
-To illustrate, consider the interaction between clustering and $t$-SNE.
-We do not perform clustering on the $t$-SNE coordinates, but rather, we cluster on the first 10-50 PCs (Chapter \@ref(clustering)) and then visualize the cluster identities on $t$-SNE plots like that in Figure \@ref(fig:tsne-brain).
-This ensures that clustering makes use of the information that was lost during compression into two dimensions for visualization.
-The plot can then be used for a diagnostic inspection of the clustering output, e.g., to check which clusters are close neighbors or whether a cluster can be split into further subclusters; this follows the aforementioned theme of focusing on local structure.
-
-From a naive perspective, using the $t$-SNE coordinates directly for clustering is tempting as it ensures that any results are immediately consistent with the visualization.
-Given that clustering is rather arbitrary anyway, there is nothing inherently wrong with this strategy - in fact, it can be treated as a rather circuitous implementation of graph-based clustering (Section \@ref(clustering-graph)).
-However, the enforced consistency can actually be considered a disservice as it masks the ambiguity of the conclusions, either due to the loss of information from dimensionality reduction or the uncertainty of the clustering.
-Rather than being errors, major discrepancies can instead be useful for motivating further investigation into the less obvious aspects of the dataset; conversely, the lack of discrepancies increases trust in the conclusions.
-
-## Session Info {-}
-
-```{r sessionInfo, echo=FALSE, results='asis'}
-prettySessionInfo()
-```
diff --git a/inst/book/ref.bib b/inst/book/ref.bib
deleted file mode 100644
index bfb21de..0000000
--- a/inst/book/ref.bib
+++ /dev/null
@@ -1,1743 +0,0 @@
-@Article{stegle2015computational,
- Author="Stegle, O. and Teichmann, S. A. and Marioni, J. C. ",
- Title="{{C}omputational and analytical challenges in single-cell transcriptomics}",
- Journal="Nat. Rev. Genet.",
- Year="2015",
- Volume="16",
- Number="3",
- Pages="133--145",
- Month="Mar"
-}
-
-@Article{islam2014quantitative,
- Author="Islam, S. and Zeisel, A. and Joost, S. and La Manno, G. and Zajac, P. and Kasper, M. and Lonnerberg, P. and Linnarsson, S. ",
- Title="{{Q}uantitative single-cell {R}{N}{A}-seq with unique molecular identifiers}",
- Journal="Nat. Methods",
- Year="2014",
- Volume="11",
- Number="2",
- Pages="163--166",
- Month="Feb"
-}
-
-@Article{huber2015orchestrating,
- Author="Huber, W. and Carey, V. J. and Gentleman, R. and Anders, S. and Carlson, M. and Carvalho, B. S. and Bravo, H. C. and Davis, S. and Gatto, L. and Girke, T. and Gottardo, R. and Hahne, F. and Hansen, K. D. and Irizarry, R. A. and Lawrence, M. and Love, M. I. and MacDonald, J. and Obenchain, V. and Oles, A. K. and Pag\`es, H. and Reyes, A. and Shannon, P. and Smyth, G. K. and Tenenbaum, D. and Waldron, L. and Morgan, M. ",
- Title="{{O}rchestrating high-throughput genomic analysis with {B}ioconductor}",
- Journal="Nat. Methods",
- Year="2015",
- Volume="12",
- Number="2",
- Pages="115--121",
- Month="Feb"
-}
-
-@Article{brennecke2013accounting,
- Author="Brennecke, P. and Anders, S. and Kim, J. K. and Ko\l{}odziejczyk, A. A. and Zhang, X. and Proserpio, V. and Baying, B. and Benes, V. and Teichmann, S. A. and Marioni, J. C. and Heisler, M. G. ",
- Title="{{A}ccounting for technical noise in single-cell {R}{N}{A}-seq experiments}",
- Journal="Nat. Methods",
- Year="2013",
- Volume="10",
- Number="11",
- Pages="1093--1095",
- Month="Nov"
-}
-
-@Article{marinov2014singlecell,
- Author="Marinov, G. K. and Williams, B. A. and McCue, K. and Schroth, G. P. and Gertz, J. and Myers, R. M. and Wold, B. J. ",
- Title="{{F}rom single-cell to cell-pool transcriptomes: stochasticity in gene expression and {R}{N}{A} splicing}",
- Journal="Genome Res.",
- Year="2014",
- Volume="24",
- Number="3",
- Pages="496--510",
- Month="Mar"
-}
-
-@Article{vallejos2015basics,
- Author="Vallejos, C. A. and Marioni, J. C. and Richardson, S. ",
- Title="{{B}{A}{S}i{C}{S}: Bayesian analysis of single-cell sequencing data}",
- Journal="PLoS Comput. Biol.",
- Year="2015",
- Volume="11",
- Number="6",
- Pages="e1004333",
- Month="Jun"
-}
-
-@Article{trapnell2014dynamics,
- Author="Trapnell, C. and Cacchiarelli, D. and Grimsby, J. and Pokharel, P. and Li, S. and Morse, M. and Lennon, N. J. and Livak, K. J. and Mikkelsen, T. S. and Rinn, J. L. ",
- Title="{{T}he dynamics and regulators of cell fate decisions are revealed by pseudotemporal ordering of single cells}",
- Journal="Nat. Biotechnol.",
- Year="2014",
- Volume="32",
- Number="4",
- Pages="381--386",
- Month="Apr"
-}
-
-@Article{kolod2015singlecell,
- Author="Ko\l{}odziejczyk, A. A. and Kim, J. K. and Tsang, J. C. and Ilicic, T. and Henriksson, J. and Natarajan, K. N. and Tuck, A. C. and Gao, X. and Buhler, M. and Liu, P. and Marioni, J. C. and Teichmann, S. A. ",
- Title="{{S}ingle cell {R}{N}{A}-sequencing of pluripotent states unlocks modular transcriptional variation}",
- Journal="Cell Stem Cell",
- Year="2015",
- Volume="17",
- Number="4",
- Pages="471--485",
- Month="Oct"
-}
-
-@Article{anders2010differential,
- Author="Anders, S. and Huber, W. ",
- Title="{{D}ifferential expression analysis for sequence count data}",
- Journal="Genome Biol.",
- Year="2010",
- Volume="11",
- Number="10",
- Pages="R106"
-}
-
-@Article{robinson2010scaling,
- Author="Robinson, M. D. and Oshlack, A. ",
- Title="{{A} scaling normalization method for differential expression analysis of {R}{N}{A}-seq data}",
- Journal="Genome Biol.",
- Year="2010",
- Volume="11",
- Number="3",
- Pages="R25"
-}
-
-@Article{islam2011characterization,
- Author="Islam, S. and Kjallquist, U. and Moliner, A. and Zajac, P. and Fan, J. B. and Lonnerberg, P. and Linnarsson, S. ",
- Title="{{C}haracterization of the single-cell transcriptional landscape by highly multiplex {R}{N}{A}-seq}",
- Journal="Genome Res.",
- Year="2011",
- Volume="21",
- Number="7",
- Pages="1160--1167",
- Month="Jul"
-}
-
-@Article{wilson2015combined,
- Author="Wilson, N. K. and Kent, D. G. and Buettner, F. and Shehata, M. and Macaulay, I. C. and Calero-Nieto, F. J. and Sanchez Castillo, M. and Oedekoven, C. A. and Diamanti, E. and Schulte, R. and Ponting, C. P. and Voet, T. and Caldas, C. and Stingl, J. and Green, A. R. and Theis, F. J. and Gottgens, B. ",
- Title="{{C}ombined single-cell functional and gene expression analysis resolves heterogeneity within stem cell populations}",
- Journal="Cell Stem Cell",
- Year="2015",
- Volume="16",
- Number="6",
- Pages="712--724",
- Month="Jun"
-}
-
-@Article{picelli2014fulllength,
- Author="Picelli, S. and Faridani, O. R. and Bjorklund, A. K. and Winberg, G. and Sagasser, S. and Sandberg, R. ",
- Title="{{F}ull-length {R}{N}{A}-seq from single cells using {S}mart-seq2}",
- Journal="Nat Protoc",
- Year="2014",
- Volume="9",
- Number="1",
- Pages="171--181",
- Month="Jan"
-}
-
-@Article{pollen2014lowcoverage,
- Author="Pollen, A. A. and Nowakowski, T. J. and Shuga, J. and Wang, X. and Leyrat, A. A. and Lui, J. H. and Li, N. and Szpankowski, L. and Fowler, B. and Chen, P. and Ramalingam, N. and Sun, G. and Thu, M. and Norris, M. and Lebofsky, R. and Toppani, D. and Kemp, D. W. and Wong, M. and Clerkson, B. and Jones, B. N. and Wu, S. and Knutsson, L. and Alvarado, B. and Wang, J. and Weaver, L. S. and May, A. P. and Jones, R. C. and Unger, M. A. and Kriegstein, A. R. and West, J. A. ",
- Title="{{L}ow-coverage single-cell m{R}{N}{A} sequencing reveals cellular heterogeneity and activated signaling pathways in developing cerebral cortex}",
- Journal="Nat. Biotechnol.",
- Year="2014",
- Volume="32",
- Number="10",
- Pages="1053--1058",
- Month="Oct"
-}
-
-@Article{love2015rnaseq,
- Author="Love, M. I. and Anders, S. and Kim, V. and Huber, W. ",
- Title="{{R}{N}{A}-{S}eq workflow: gene-level exploratory analysis and differential expression}",
- Journal="F1000Res",
- Year="2015",
- Volume="4",
- Pages="1070"
-}
-
-@Article{liao2014featurecounts,
- Author="Liao, Y. and Smyth, G. K. and Shi, W. ",
- Title="{feature{C}ounts: an efficient general purpose program for assigning sequence reads to genomic features}",
- Journal="Bioinformatics",
- Year="2014",
- Volume="30",
- Number="7",
- Pages="923--930",
- Month="Apr"
-}
-
-@Article{liao2013subread,
- Author="Liao, Y. and Smyth, G. K. and Shi, W. ",
- Title="{{T}he {S}ubread aligner: fast, accurate and scalable read mapping by seed-and-vote}",
- Journal="Nucleic Acids Res.",
- Year="2013",
- Volume="41",
- Number="10",
- Pages="e108",
- Month="May"
-}
-
-@Article{law2014voom,
- Author="Law, C. W. and Chen, Y. and Shi, W. and Smyth, G. K. ",
- Title="{voom: {P}recision weights unlock linear model analysis tools for {R}{N}{A}-seq read counts}",
- Journal="Genome Biol.",
- Year="2014",
- Volume="15",
- Number="2",
- Pages="R29"
-}
-
-@Article{love2014moderated,
- Author="Love, M. I. and Huber, W. and Anders, S. ",
- Title="{{M}oderated estimation of fold change and dispersion for {R}{N}{A}-seq data with {D}{E}{S}eq2}",
- Journal="Genome Biol.",
- Year="2014",
- Volume="15",
- Number="12",
- Pages="550"
-}
-
-@article{van2008visualizing,
- title={Visualizing data using t-{S}{N}{E}},
- author={Van der Maaten, L. and Hinton, G.},
- journal={J. Mach. Learn. Res.},
- volume={9},
- number={2579-2605},
- pages={85},
- year={2008}
-}
-
-@Article{scialdone2015computational,
- Author="Scialdone, A. and Natarajan, K. N. and Saraiva, L. R. and Proserpio, V. and Teichmann, S. A. and Stegle, O. and Marioni, J. C. and Buettner, F. ",
- Title="{{C}omputational assignment of cell-cycle stage from single-cell transcriptome data}",
- Journal="Methods",
- Year="2015",
- Volume="85",
- Pages="54--61",
- Month="Sep"
-}
-
-@Article{kim2015characterizing,
- Author="Kim, J. K. and Ko\l{}odziejczyk, A. A. and Illicic, T. and Teichmann, S. A. and Marioni, J. C. ",
- Title="{{C}haracterizing noise structure in single-cell {R}{N}{A}-seq distinguishes genuine from technical stochastic allelic expression}",
- Journal="Nat. Commun.",
- Year="2015",
- Volume="6",
- Pages="8687"
-}
-
-@Article{mccarthy2012differential,
- Author="McCarthy, D. J. and Chen, Y. and Smyth, G. K. ",
- Title="{{D}ifferential expression analysis of multifactor {R}{N}{A}-{S}eq experiments with respect to biological variation}",
- Journal="Nucleic Acids Res.",
- Year="2012",
- Volume="40",
- Number="10",
- Pages="4288--4297",
- Month="May"
-}
-
-@Article{angel1991role,
- Author="Angel, P. and Karin, M. ",
- Title="{{T}he role of {J}un, {F}os and the {A}{P}-1 complex in cell-proliferation and transformation}",
- Journal="Biochim. Biophys. Acta",
- Year="1991",
- Volume="1072",
- Number="2-3",
- Pages="129--157",
- Month="Dec"
-}
-
-@Article{langfelder2008defining,
- Author="Langfelder, P. and Zhang, B. and Horvath, S. ",
- Title="{{D}efining clusters from a hierarchical cluster tree: the {D}ynamic {T}ree {C}ut package for {R}}",
- Journal="Bioinformatics",
- Year="2008",
- Volume="24",
- Number="5",
- Pages="719--720",
- Month="Mar"
-}
-
-@Article{ilicic2016classification,
- Author="Ilicic, T. and Kim, J. K. and Ko\l{}odziejczyk, A. A. and Bagger, F. O. and McCarthy, D. J. and Marioni, J. C. and Teichmann, S. A. ",
- Title="{{C}lassification of low quality cells from single-cell {R}{N}{A}-seq data}",
- Journal="Genome Biol.",
- Year="2016",
- Volume="17",
- Number="1",
- Pages="29"
-}
-
-@Article{heng2008immunological,
- Author="Heng, T. S. and Painter, M. W. and Elpek, K. and Lukacs-Kornek, V. and Mauermann, N. and Turley, S. J. and Koller, D. and Kim, F. S. and Wagers, A. J. and Asinovski, N. and Davis, S. and Fassett, M. and Feuerer, M. and Gray, D. H. and Haxhinasto, S. and Hill, J. A. and Hyatt, G. and Laplace, C. and Leatherbee, K. and Mathis, D. and Benoist, C. and Jianu, R. and Laidlaw, D. H. and Best, J. A. and Knell, J. and Goldrath, A. W. and Jarjoura, J. and Sun, J. C. and Zhu, Y. and Lanier, L. L. and Ergun, A. and Li, Z. and Collins, J. J. and Shinton, S. A. and Hardy, R. R. and Friedline, R. and Sylvia, K. and Kang, J. ",
- Title="{{T}he {I}mmunological {G}enome {P}roject: networks of gene expression in immune cells}",
- Journal="Nat. Immunol.",
- Year="2008",
- Volume="9",
- Number="10",
- Pages="1091--1094",
- Month="Oct"
-}
-
-@Article{zeisel2015brain,
- Author="Zeisel, A. and Munoz-Manchado, A. B. and Codeluppi, S. and Lonnerberg, P. and La Manno, G. and Jureus, A. and Marques, S. and Munguba, H. and He, L. and Betsholtz, C. and Rolny, C. and Castelo-Branco, G. and Hjerling-Leffler, J. and Linnarsson, S. ",
- Title="{{B}rain structure. {C}ell types in the mouse cortex and hippocampus revealed by single-cell {R}{N}{A}-seq}",
- Journal="Science",
- Year="2015",
- Volume="347",
- Number="6226",
- Pages="1138--1142",
- Month="Mar"
-}
-
-@Article{robinson2010edgeR,
- Author="Robinson, M. D. and McCarthy, D. J. and Smyth, G. K. ",
- Title="{edge{R}: a {B}ioconductor package for differential expression analysis of digital gene expression data}",
- Journal="Bioinformatics",
- Year="2010",
- Volume="26",
- Number="1",
- Pages="139--140",
- Month="Jan"
-}
-
-@Article{ritchie2015limma,
- Author="Ritchie, M. E. and Phipson, B. and Wu, D. and Hu, Y. and Law, C. W. and Shi, W. and Smyth, G. K. ",
- Title="{limma powers differential expression analyses for {R}{N}{A}-sequencing and microarray studies}",
- Journal="Nucleic Acids Res.",
- Year="2015",
- Volume="43",
- Number="7",
- Pages="e47",
- Month="Apr"
-}
-
-@Article{julia2015sincell,
- Author="Julia, M. and Telenti, A. and Rausell, A. ",
- Title="{{S}incell: an {R}/{B}ioconductor package for statistical assessment of cell-state hierarchies from single-cell {R}{N}{A}-seq}",
- Journal="Bioinformatics",
- Year="2015",
- Volume="31",
- Number="20",
- Pages="3380--3382",
- Month="Oct"
-}
-
-@Article{buettner2015computational,
- Author="Buettner, F. and Natarajan, K. N. and Casale, F. P. and Proserpio, V. and Scialdone, A. and Theis, F. J. and Teichmann, S. A. and Marioni, J. C. and Stegle, O. ",
- Title="{{C}omputational analysis of cell-to-cell heterogeneity in single-cell {R}{N}{A}-sequencing data reveals hidden subpopulations of cells}",
- Journal="Nat. Biotechnol.",
- Year="2015",
- Volume="33",
- Number="2",
- Pages="155--160",
- Month="Feb"
-}
-
-@Article{mahata2014singlecell,
- Author="Mahata, B. and Zhang, X. and Ko\l{}odziejczyk, A. A. and Proserpio, V. and Haim-Vilmovsky, L. and Taylor, A. E. and Hebenstreit, D. and Dingler, F. A. and Moignard, V. and Gottgens, B. and Arlt, W. and McKenzie, A. N. and Teichmann, S. A. ",
- Title="{{S}ingle-cell {R}{N}{A} sequencing reveals {T} helper cells synthesizing steroids de novo to contribute to immune homeostasis}",
- Journal="Cell Rep.",
- Year="2014",
- Volume="7",
- Number="4",
- Pages="1130--1142",
- Month="May"
-}
-
-@Article{phipson2010permutation,
- Author="Phipson, B. and Smyth, G. K. ",
- Title="{{P}ermutation {P}-values should never be zero: calculating exact {P}-values when permutations are randomly drawn}",
- Journal="Stat. Appl. Genet. Mol. Biol.",
- Year="2010",
- Volume="9",
- Pages="Article 39"
-}
-
-@incollection{chen2014differential,
- title="{Differential expression analysis of complex RNA-seq experiments using edgeR}",
- author={Chen, Y. and Lun, A. T. L. and Smyth, G. K.},
- booktitle={Statistical analysis of next generation sequencing data},
- editor={Datta, S. and Nettleton, D.},
- pages={51--74},
- year={2014},
- publisher={Springer}
-}
-
-@article{hicks2015widespread,
- author = {Hicks, S. C. and Teng, M. and Irizarry, R. A.},
- title = "{On the widespread and critical impact of systematic bias and batch effects in single-cell RNA-Seq data}",
- year = {2015},
- doi = {10.1101/025528},
- publisher = {Cold Spring Harbor Labs Journals},
- URL = {http://biorxiv.org/content/early/2015/08/25/025528},
- eprint = {http://biorxiv.org/content/early/2015/08/25/025528.full.pdf},
- journal = {bioRxiv}
-}
-
-@Article{leng2015oscope,
- Author="Leng, N. and Chu, L. F. and Barry, C. and Li, Y. and Choi, J. and Li, X. and Jiang, P. and Stewart, R. M. and Thomson, J. A. and Kendziorski, C. ",
- Title="{{O}scope identifies oscillatory genes in unsynchronized single-cell {R}{N}{A}-seq experiments}",
- Journal="Nat. Methods",
- Year="2015",
- Volume="12",
- Number="10",
- Pages="947--950",
- Month="Oct"
-}
-
-@Article{fan2016characterizing,
- Author="Fan, J. and Salathia, N. and Liu, R. and Kaeser, G. E. and Yung, Y. C. and Herman, J. L. and Kaper, F. and Fan, J. B. and Zhang, K. and Chun, J. and Kharchenko, P. V. ",
- Title="{{C}haracterizing transcriptional heterogeneity through pathway and gene set overdispersion analysis}",
- Journal="Nat. Methods",
- Year="2016",
- Volume="13",
- Number="3",
- Pages="241--244",
- Month="Mar"
-}
-
-@Article{lund2012detecting,
- Author="Lund, S. P. and Nettleton, D. and McCarthy, D. J. and Smyth, G. K. ",
- Title="{{D}etecting differential expression in {R}{N}{A}-sequence data using quasi-likelihood with shrunken dispersion estimates}",
- Journal="Stat Appl Genet Mol Biol",
- Year="2012",
- Volume="11",
- Number="5",
- Pages=" "
-}
-
-@Article{angerer2016destiny,
- Author="Angerer, P. and Haghverdi, L. and Buttner, M. and Theis, F. J. and Marr, C. and Buettner, F. ",
- Title="{destiny: diffusion maps for large-scale single-cell data in {R}}",
- Journal="Bioinformatics",
- Year="2016",
- Volume="32",
- Number="8",
- Pages="1241--1243",
-}
-
-@Article{lun2016pooling,
- Author="Lun, A. T. and Bach, K. and Marioni, J. C. ",
- Title="{{P}ooling across cells to normalize single-cell {R}{N}{A} sequencing data with many zero counts}",
- Journal="Genome Biol.",
- Year="2016",
- Volume="17",
- Pages="75",
- Month="Apr"
-}
-
-@Article{zhu2006gata3,
- Author="Zhu, J. and Yamane, H. and Cote-Sierra, J. and Guo, L. and Paul, W. E.",
- Title="{{G}{A}{T}{A}-3 promotes {T}h2 responses through three different mechanisms: induction of {T}h2 cytokine production, selective growth of {T}h2 cells and inhibition of {T}h1 cell-specific factors}",
- Journal="Cell Res.",
- Year="2006",
- Volume="16",
- Number="1",
- Pages="3--10",
- Month="Jan"
-}
-
-@article{berger1996bioequivalence,
- author = "Berger, R. L. and Hsu, J. C.",
- doi = "10.1214/ss/1032280304",
- fjournal = "Statistical Science",
- journal = "Statist. Sci.",
- month = "11",
- number = "4",
- pages = "283--319",
- publisher = "The Institute of Mathematical Statistics",
- title = "Bioequivalence trials, intersection-union tests and equivalence confidence sets",
- url = "http://dx.doi.org/10.1214/ss/1032280304",
- volume = "11",
- year = "1996"
-}
-
-@Article{santos2015cyclebase,
- Author="Santos, A. and Wernersson, R. and Jensen, L. J. ",
- Title="{{C}yclebase 3.0: a multi-organism database on cell-cycle regulation and phenotypes}",
- Journal="Nucleic Acids Res.",
- Year="2015",
- Volume="43",
- Number="Database issue",
- Pages="D1140--1144",
- Month="Jan"
-}
-
-@article {tung2017batch,
- Author="Tung, P. Y. and Blischak, J. D. and Hsiao, C. J. and Knowles, D. A. and Burnett, J. E. and Pritchard, J. K. and Gilad, Y. ",
- Title="{{B}atch effects and the effective design of single-cell gene expression studies}",
- Journal="Sci. Rep.",
- Year="2017",
- Volume="7",
- Pages="39921",
- Month="Jan"
-}
-
-@Article{bray2016near,
- Author="Bray, N. L. and Pimentel, H. and Melsted, P. and Pachter, L. ",
- Title="{{N}ear-optimal probabilistic {R}{N}{A}-seq quantification}",
- Journal="Nat. Biotechnol.",
- Year="2016",
- Volume="34",
- Number="5",
- Pages="525--527",
- Month="May"
-}
-
-@article{patro2015accurate,
- author = {Patro, Rob and Duggal, Geet and Kingsford, Carl},
- title = {Accurate, fast, and model-aware transcript expression quantification with Salmon},
- year = {2015},
- doi = {10.1101/021592},
- publisher = {Cold Spring Harbor Labs Journals},
- journal = {bioRxiv}
-}
-
-@Article{ji2016tscan,
- Author="Ji, Z. and Ji, H. ",
- Title="{{T}{S}{C}{A}{N}: {P}seudo-time reconstruction and evaluation in single-cell {R}{N}{A}-seq analysis}",
- Journal="Nucleic Acids Res.",
- Year="2016",
- Volume="44",
- Number="13",
- Pages="e117",
- Month="Jul"
-}
-
-@Article{klein2015droplet,
- Author="Klein, A. M. and Mazutis, L. and Akartuna, I. and Tallapragada, N. and Veres, A. and Li, V. and Peshkin, L. and Weitz, D. A. and Kirschner, M. W. ",
- Title="{{D}roplet barcoding for single-cell transcriptomics applied to embryonic stem cells}",
- Journal="Cell",
- Year="2015",
- Volume="161",
- Number="5",
- Pages="1187--1201",
- Month="May"
-}
-
-@Article{macosko2015highly,
- Author="Macosko, E. Z. and Basu, A. and Satija, R. and Nemesh, J. and Shekhar, K. and Goldman, M. and Tirosh, I. and Bialas, A. R. and Kamitaki, N. and Martersteck, E. M. and Trombetta, J. J. and Weitz, D. A. and Sanes, J. R. and Shalek, A. K. and Regev, A. and McCarroll, S. A. ",
- Title="{{H}ighly parallel genome-wide expression profiling of individual cells using nanoliter droplets}",
- Journal="Cell",
- Year="2015",
- Volume="161",
- Number="5",
- Pages="1202--1214",
- Month="May"
-}
-
-@Article{bourgon2010independent,
- Author="Bourgon, R. and Gentleman, R. and Huber, W. ",
- Title="{{I}ndependent filtering increases detection power for high-throughput experiments}",
- Journal="Proc. Natl. Acad. Sci. U.S.A.",
- Year="2010",
- Volume="107",
- Number="21",
- Pages="9546--9551",
- Month="May"
-}
-
-@Article{finak2015mast,
- Author="Finak, G. and McDavid, A. and Yajima, M. and Deng, J. and Gersuk, V. and Shalek, A. K. and Slichter, C. K. and Miller, H. W. and McElrath, M. J. and Prlic, M. and Linsley, P. S. and Gottardo, R. ",
- Title="{{M}{A}{S}{T}: a flexible statistical framework for assessing transcriptional changes and characterizing heterogeneity in single-cell {R}{N}{A} sequencing data}",
- Journal="Genome Biol.",
- Year="2015",
- Volume="16",
- Pages="278"
-}
-
-@Article{conboy2007cell,
- Author="Conboy, C. M. and Spyrou, C. and Thorne, N. P. and Wade, E. J. and Barbosa-Morais, N. L. and Wilson, M. D. and Bhattacharjee, A. and Young, R. A. and Tavare, S. and Lees, J. A. and Odom, D. T. ",
- Title="{{C}ell cycle genes are the evolutionarily conserved targets of the {E}2{F}4 transcription factor}",
- Journal="PLoS ONE",
- Year="2007",
- Volume="2",
- Number="10",
- Pages="e1061"
-}
-
-@Article{bertoli2013control,
- Author="Bertoli, C. and Skotheim, J. M. and de Bruin, R. A. ",
- Title="{{C}ontrol of cell cycle transcription during {G}1 and {S} phases}",
- Journal="Nat. Rev. Mol. Cell Biol.",
- Year="2013",
- Volume="14",
- Number="8",
- Pages="518--528",
- Month="Aug"
-}
-
-@Article{chen2016from,
- Author="Chen, Y. and Lun, A. T. and Smyth, G. K. ",
- Title="{{F}rom reads to genes to pathways: differential expression analysis of {R}{N}{A}-{S}eq experiments using {R}subread and the edge{R} quasi-likelihood pipeline}",
- Journal="F1000Res",
- Year="2016",
- Volume="5",
- Pages="1438"
-}
-
-@Article{coller2006new,
- Author="Coller, H. A. and Sang, L. and Roberts, J. M. ",
- Title="{{A} new description of cellular quiescence}",
- Journal="PLoS Biol.",
- Year="2006",
- Volume="4",
- Number="3",
- Pages="e83",
- Month="Mar"
-}
-
-@Article{mccarthy2017scater,
- Author="McCarthy, D. J. and Campbell, K. R. and Lun, A. T. and Wills, Q. F. ",
- Title="{{S}cater: pre-processing, quality control, normalization and visualization of single-cell {R}{N}{A}-seq data in {R}}",
- Journal="Bioinformatics",
- Year="2017",
- Volume="33",
- Number="8",
- Pages="1179--1186",
- Month="Apr"
-}
-
-@article{simes1986improved,
- author = {Simes, R. J.},
- title = {An improved {B}onferroni procedure for multiple tests of significance},
- volume = {73},
- number = {3},
- pages = {751--754},
- year = {1986},
- journal = {Biometrika}
-}
-
-@Article{xu2015identification,
- Author="Xu, C. and Su, Z. ",
- Title="{{I}dentification of cell types from single-cell transcriptomes using a novel clustering method}",
- Journal="Bioinformatics",
- Year="2015",
- Volume="31",
- Number="12",
- Pages="1974--1980",
- Month="Jun"
-}
-
-@ARTICLE{lun2016step,
- title = "A step-by-step workflow for low-level analysis of single-cell
- {RNA-seq} data",
- author = "Lun, A. T. L. and McCarthy, D. J. and Marioni, J. C.",
- journal = "F1000Res.",
- volume = 5,
- month = aug,
- year = 2016
-}
-
-@Article{lun2017assessing,
- Author="Lun, A. T. L. and Calero-Nieto, F. J. and Haim-Vilmovsky, L. and Gottgens, B. and Marioni, J. C. ",
- Title="{{A}ssessing the reliability of spike-in normalization for analyses of single-cell {R}{N}{A} sequencing data}",
- Journal="Genome Res.",
- Year="2017",
- Volume="27",
- Number="11",
- Pages="1795--1806",
- Month="Nov"
-}
-
-@Article{lun2017overcoming,
- Author="Lun, A. T. L. and Marioni, J. C. ",
- Title="{{O}vercoming confounding plate effects in differential expression analyses of single-cell {R}{N}{A}-seq data}",
- Journal="Biostatistics",
- Year="2017",
- Volume="18",
- Number="3",
- Pages="451--464",
- Month="Jul"
-}
-
-@Article{yang2016comparative,
- Author="Yang, Z. and Algesheimer, R. and Tessone, C. J. ",
- Title="{A comparative analysis of community detection algorithms on artificial networks}",
- Journal="Sci. Rep.",
- Year="2016",
- Volume="6",
- Pages="30750",
- Month="Aug"
-}
-
-@Article{smith2017umitools,
- Author="Smith, T. and Heger, A. and Sudbery, I. ",
- Title="{{U}{M}{I}-tools: modeling sequencing errors in {U}nique {M}olecular {I}dentifiers to improve quantification accuracy}",
- Journal="Genome Res.",
- Year="2017",
- Volume="27",
- Number="3",
- Pages="491--499",
- Month="Mar"
-}
-
-@Article{tasic2016adult,
- Author="Tasic, B. and Menon, V. and Nguyen, T. N. and Kim, T. K. and Jarsky, T. and Yao, Z. and Levi, B. and Gray, L. T. and Sorensen, S. A. and Dolbeare, T. and Bertagnolli, D. and Goldy, J. and Shapovalova, N. and Parry, S. and Lee, C. and Smith, K. and Bernard, A. and Madisen, L. and Sunkin, S. M. and Hawrylycz, M. and Koch, C. and Zeng, H. ",
- Title="{{A}dult mouse cortical cell taxonomy revealed by single cell transcriptomics}",
- Journal="Nat. Neurosci.",
- Year="2016",
- Volume="19",
- Number="2",
- Pages="335--346",
- Month="Feb"
-}
-
-@Article{wajapeyee2010senescence,
- Author="Wajapeyee, N. and Wang, S. Z. and Serra, R. W. and Solomon, P. D. and Nagarajan, A. and Zhu, X. and Green, M. R. ",
- Title="{{S}enescence induction in human fibroblasts and hematopoietic progenitors by leukemogenic fusion proteins}",
- Journal="Blood",
- Year="2010",
- Volume="115",
- Number="24",
- Pages="5057--5060",
- Month="Jun"
-}
-
-@Article{zeng2012largescale,
- Author="Zeng, H. and Shen, E. H. and Hohmann, J. G. and Oh, S. W. and Bernard, A. and Royall, J. J. and Glattfelder, K. J. and Sunkin, S. M. and Morris, J. A. and Guillozet-Bongaarts, A. L. and Smith, K. A. and Ebbert, A. J. and Swanson, B. and Kuan, L. and Page, D. T. and Overly, C. C. and Lein, E. S. and Hawrylycz, M. J. and Hof, P. R. and Hyde, T. M. and Kleinman, J. E. and Jones, A. R. ",
- Title="{{L}arge-scale cellular-resolution gene profiling in human neocortex reveals species-specific molecular signatures}",
- Journal="Cell",
- Year="2012",
- Volume="149",
- Number="2",
- Pages="483--496",
- Month="Apr"
-}
-
-@Article{grun2016denovo,
- Author="Grun, D. and Muraro, M. J. and Boisset, J. C. and Wiebrands, K. and Lyubimova, A. and Dharmadhikari, G. and van den Born, M. and van Es, J. and Jansen, E. and Clevers, H. and de Koning, E. J. P. and van Oudenaarden, A. ",
- Title="{{D}e {N}ovo {P}rediction of {S}tem {C}ell {I}dentity using {S}ingle-{C}ell {T}ranscriptome {D}ata}",
- Journal="Cell Stem Cell",
- Year="2016",
- Volume="19",
- Number="2",
- Pages="266--277",
- Month="08"
-}
-
-@article{haghverdi2018batch,
- Author="Haghverdi, L. and Lun, A. T. L. and Morgan, M. D. and Marioni, J. C. ",
- Title="{{B}atch effects in single-cell {R}{N}{A}-sequencing data are corrected by matching mutual nearest neighbors}",
- Journal="Nat. Biotechnol.",
- Year="2018",
- Volume="36",
- Number="5",
- Pages="421--427"
-}
-
-@Article{muraro2016singlecell,
- Author="Muraro, M. J. and Dharmadhikari, G. and Grun, D. and Groen, N. and Dielen, T. and Jansen, E. and van Gurp, L. and Engelse, M. A. and Carlotti, F. and de Koning, E. J. and van Oudenaarden, A. ",
- Title="{{A} {S}ingle-{C}ell {T}ranscriptome {A}tlas of the {H}uman {P}ancreas}",
- Journal="Cell Syst",
- Year="2016",
- Volume="3",
- Number="4",
- Pages="385--394",
- Month="Oct"
-}
-
-@Article{segerstolpe2016singlecell,
- Author="Segerstolpe, A. and Palasantza, A. and Eliasson, P. and Andersson, E. M. and Andreasson, A. C. and Sun, X. and Picelli, S. and Sabirsh, A. and Clausen, M. and Bjursell, M. K. and Smith, D. M. and Kasper, M. and Ammala, C. and Sandberg, R. ",
- Title="{{S}ingle-{C}ell {T}ranscriptome {P}rofiling of {H}uman {P}ancreatic {I}slets in {H}ealth and {T}ype 2 {D}iabetes}",
- Journal="Cell Metab.",
- Year="2016",
- Volume="24",
- Number="4",
- Pages="593--607",
- Month="Oct"
-}
-
-@Article{horn1965rationale,
- author="Horn, J. L.",
- title="A rationale and test for the number of factors in factor analysis",
- journal="Psychometrika",
- year="1965",
- month="Jun",
- day="01",
- volume="30",
- number="2",
- pages="179--185",
- issn="1860-0980",
- doi="10.1007/BF02289447",
- url="https://doi.org/10.1007/BF02289447"
-}
-
-@Article{buja1992remarks,
- Author="Buja, A. and Eyuboglu, N. ",
- Title="{{R}emarks on {P}arallel {A}nalysis}",
- Journal="Multivariate Behav Res",
- Year="1992",
- Volume="27",
- Number="4",
- Pages="509--540",
- Month="Oct"
-}
-
-@Article{zheng2017massively,
- Author="Zheng, G. X. and Terry, J. M. and Belgrader, P. and Ryvkin, P. and Bent, Z. W. and Wilson, R. and Ziraldo, S. B. and Wheeler, T. D. and McDermott, G. P. and Zhu, J. and Gregory, M. T. and Shuga, J. and Montesclaros, L. and Underwood, J. G. and Masquelier, D. A. and Nishimura, S. Y. and Schnall-Levin, M. and Wyatt, P. W. and Hindson, C. M. and Bharadwaj, R. and Wong, A. and Ness, K. D. and Beppu, L. W. and Deeg, H. J. and McFarland, C. and Loeb, K. R. and Valente, W. J. and Ericson, N. G. and Stevens, E. A. and Radich, J. P. and Mikkelsen, T. S. and Hindson, B. J. and Bielas, J. H. ",
- Title="{{M}assively parallel digital transcriptional profiling of single cells}",
- Journal="Nat Commun",
- Year="2017",
- Volume="8",
- Pages="14049",
- Month="Jan"
-}
-
-@article{lun2018distinguishing,
- author = {Lun, A. and Riesenfeld, S. and Andrews, T. and Dao, T. P. and Gomes, T. and {{participants in the 1\textsuperscript{st} Human Cell Atlas Jamboree}} and Marioni, J.},
- Title="{{E}mpty{D}rops: distinguishing cells from empty droplets in droplet-based single-cell {R}{N}{A} sequencing data}",
- Journal="Genome Biol.",
- Year="2019",
- Volume="20",
- Number="1",
- Pages="63",
- Month="03"
-}
-
-@article {wolf2017graph,
- author = {Wolf, F. Alexander and Hamey, Fiona and Plass, Mireya and Solana, Jordi and Dahlin, Joakim S. and Gottgens, Berthold and Rajewsky, Nikolaus and Simon, Lukas and Theis, Fabian J.},
- title = {Graph abstraction reconciles clustering with trajectory inference through a topology preserving map of single cells},
- year = {2017},
- doi = {10.1101/208819},
- publisher = {Cold Spring Harbor Laboratory},
- URL = {https://www.biorxiv.org/content/early/2017/10/25/208819},
- eprint = {https://www.biorxiv.org/content/early/2017/10/25/208819.full.pdf},
- journal = {bioRxiv}
-}
-
-@Article{lawlor2017singlecell,
- Author="Lawlor, N. and George, J. and Bolisetty, M. and Kursawe, R. and Sun, L. and Sivakamasundari, V. and Kycia, I. and Robson, P. and Stitzel, M. L. ",
- Title="{{S}ingle-cell transcriptomes identify human islet cell signatures and reveal cell-type-specific expression changes in type 2 diabetes}",
- Journal="Genome Res.",
- Year="2017",
- Volume="27",
- Number="2",
- Pages="208--222",
- Month="02"
-}
-
-@Article{guimaraes2016patterns,
- Author="Guimaraes, J. C. and Zavolan, M. ",
- Title="{{P}atterns of ribosomal protein expression specify normal and malignant human cells}",
- Journal="Genome Biol.",
- Year="2016",
- Volume="17",
- Number="1",
- Pages="236",
- Month="11"
-}
-
-@Article{glare2002betaactin,
- Author="Glare, E. M. and Divjak, M. and Bailey, M. J. and Walters, E. H. ",
- Title="{beta-{A}ctin and {G}{A}{P}{D}{H} housekeeping gene expression in asthmatic airways is variable and not suitable for normalising m{R}{N}{A} levels}",
- Journal="Thorax",
- Year="2002",
- Volume="57",
- Number="9",
- Pages="765--770",
- Month="Sep"
-}
-
-@Article{nazari2015gapdh,
- Author="Nazari, F. and Parham, A. and Maleki, A. F. ",
- Title="{{G}{A}{P}{D}{H}, \beta-actin and \beta-microglobulin, as three common reference genes, are not reliable for gene expression studies in equine adipose- and marrow-derived mesenchymal stem cells}",
- Journal="J Anim Sci Technol",
- Year="2015",
- Volume="57",
- Pages="18"
-}
-
-@Article{soneson2018bias,
- Author="Soneson, C. and Robinson, M. D. ",
- Title="{{B}ias, robustness and scalability in single-cell differential expression analysis}",
- Journal="Nat. Methods",
- Year="2018",
- Volume="15",
- Number="4",
- Pages="255--261",
- Month="Apr"
-}
-
-@Article{soneson2015differential,
- Author="Soneson, C. and Love, M. I. and Robinson, M. D. ",
- Title="{{D}ifferential analyses for {R}{N}{A}-seq: transcript-level estimates improve gene-level inferences}",
- Journal="F1000Res",
- Year="2015",
- Volume="4",
- Pages="1521"
-}
-
-@Article{kang2018multiplexed,
- Author="Kang, H. M. and Subramaniam, M. and Targ, S. and Nguyen, M. and Maliskova, L. and McCarthy, E. and Wan, E. and Wong, S. and Byrnes, L. and Lanata, C. M. and Gate, R. E. and Mostafavi, S. and Marson, A. and Zaitlen, N. and Criswell, L. A. and Ye, C. J. ",
- Title="{{M}ultiplexed droplet single-cell {R}{N}{A}-sequencing using natural genetic variation}",
- Journal="Nat. Biotechnol.",
- Year="2018",
- Volume="36",
- Number="1",
- Pages="89--94",
- Month="Jan"
-}
-
-@Article{dahlin2018single,
- Author="Dahlin, J. S. and Hamey, F. K. and Pijuan-Sala, B. and Shepherd, M. and Lau, W. W. Y. and Nestorowa, S. and Weinreb, C. and Wolock, S. and Hannah, R. and Diamanti, E. and Kent, D. G. and Gottgens, B. and Wilson, N. K. ",
- Title="{{A} single-cell hematopoietic landscape resolves 8 lineage trajectories and defects in {K}it mutant mice}",
- Journal="Blood",
- Year="2018",
- Volume="131",
- Number="21",
- Pages="e1-e11",
- Month="May"
-}
-
-@Article{bach2017differentiation,
- Author="Bach, K. and Pensa, S. and Grzelak, M. and Hadfield, J. and Adams, D. J. and Marioni, J. C. and Khaled, W. T. ",
- Title="{{D}ifferentiation dynamics of mammary epithelial cells revealed by single-cell {R}{N}{A} sequencing}",
- Journal="Nat Commun",
- Year="2017",
- Volume="8",
- Number="1",
- Pages="2128",
- Month="Dec"
-}
-
-@Article{lun2018beachmat,
- Author="Lun, A. T. L. and Pages, H. and Smith, M. L. ",
- Title="{beachmat: {A} {B}ioconductor {C}++ {A}{P}{I} for accessing high-throughput biological data from a variety of {R} matrix types}",
- Journal="PLoS Comput. Biol.",
- Year="2018",
- Volume="14",
- Number="5",
- Pages="e1006135",
- Month="05"
-}
-
-@Article{wang2012fast,
- author = {Wang, X.},
- title = {A fast exact k-nearest neighbors algorithm for high dimensional search using k-means clustering and triangle inequality},
- journal = {Proc Int Jt Conf Neural Netw},
- year = {2012},
- volume = {43},
- number = {6},
- pages = {2351--2358},
- month = {Feb},
-}
-
-@Article{mccarthy2009treat,
- Author="McCarthy, D. J. and Smyth, G. K. ",
- Title="{{T}esting significance relative to a fold-change threshold is a {T}{R}{E}{A}{T}}",
- Journal="Bioinformatics",
- Year="2009",
- Volume="25",
- Number="6",
- Pages="765--771",
- Month="Mar"
-}
-
-@Article{lian2018scpipe,
- Author="Tian, L. and Su, S. and Dong, X. and Amann-Zalcenstein, D. and Biben, C. and Seidi, A. and Hilton, D. J. and Naik, S. H. and Ritchie, M. E. ",
- Title="{sc{P}ipe: {A} flexible {R}/{B}ioconductor preprocessing pipeline for single-cell {R}{N}{A}-sequencing data}",
- Journal="PLoS Comput. Biol.",
- Year="2018",
- Volume="14",
- Number="8",
- Pages="e1006361",
- Month="08"
-}
-
-@Article{jaitin2014massively,
- Author="Jaitin, D. A. and Kenigsberg, E. and Keren-Shaul, H. and Elefant, N. and Paul, F. and Zaretsky, I. and Mildner, A. and Cohen, N. and Jung, S. and Tanay, A. and Amit, I. ",
- Title="{{M}assively parallel single-cell {R}{N}{A}-seq for marker-free decomposition of tissues into cell types}",
- Journal="Science",
- Year="2014",
- Volume="343",
- Number="6172",
- Pages="776--779",
- Month="Feb"
-}
-
-@Article{kang2018multiplexed,
- Author="Kang, H. M. and Subramaniam, M. and Targ, S. and Nguyen, M. and Maliskova, L. and McCarthy, E. and Wan, E. and Wong, S. and Byrnes, L. and Lanata, C. M. and Gate, R. E. and Mostafavi, S. and Marson, A. and Zaitlen, N. and Criswell, L. A. and Ye, C. J. ",
- Title="{{M}ultiplexed droplet single-cell {R}{N}{A}-sequencing using natural genetic variation}",
- Journal="Nat. Biotechnol.",
- Year="2018",
- Volume="36",
- Number="1",
- Pages="89--94",
- Month="01"
-}
-
-@Article{stoeckius2018hashing,
- Author="Stoeckius, M. and Zheng, S. and Houck-Loomis, B. and Hao, S. and Yeung, B. Z. and Mauck, W. M. and Smibert, P. and Satija, R. ",
- Title="{{C}ell {H}ashing with barcoded antibodies enables multiplexing and doublet detection for single cell genomics}",
- Journal="Genome Biol.",
- Year="2018",
- Volume="19",
- Number="1",
- Pages="224",
- Month="12"
-}
-
-@Article{davis2007geoquery,
- Author="Davis, S. and Meltzer, P. S. ",
- Title="{{G}{E}{O}query: a bridge between the {G}ene {E}xpression {O}mnibus ({G}{E}{O}) and {B}io{C}onductor}",
- Journal="Bioinformatics",
- Year="2007",
- Volume="23",
- Number="14",
- Pages="1846--1847",
- Month="Jul"
-}
-
-@Article{richard2018tcell,
- Author="Richard, A. C. and Lun, A. T. L. and Lau, W. W. Y. and Gottgens, B. and Marioni, J. C. and Griffiths, G. M. ",
- Title="{{T} cell cytolytic capacity is independent of initial stimulation strength}",
- Journal="Nat. Immunol.",
- Year="2018",
- Volume="19",
- Number="8",
- Pages="849--858",
- Month="08"
-}
-
-@article{lun2018overcoming,
- author = {Lun, A.},
- title = {Overcoming systematic errors caused by log-transformation of normalized single-cell RNA sequencing data},
- year = {2018},
- journal = {bioRxiv}
-}
-
-@Article{ji2016tscan,
- Author="Ji, Z. and Ji, H. ",
- Title="{{T}{S}{C}{A}{N}: {P}seudo-time reconstruction and evaluation in single-cell {R}{N}{A}-seq analysis}",
- Journal="Nucleic Acids Res.",
- Year="2016",
- Volume="44",
- Number="13",
- Pages="e117",
- Month="07"
-}
-
-@Article{leek2012sva,
- Author="Leek, J. T. and Johnson, W. E. and Parker, H. S. and Jaffe, A. E. and Storey, J. D. ",
- Title="{{T}he sva package for removing batch effects and other unwanted variation in high-throughput experiments}",
- Journal="Bioinformatics",
- Year="2012",
- Volume="28",
- Number="6",
- Pages="882--883",
- Month="Mar"
-}
-
-@Article{lin2019scmerge,
- Author="Lin, Y. and Ghazanfar, S. and Wang, K. Y. X. and Gagnon-Bartsch, J. A. and Lo, K. K. and Su, X. and Han, Z. G. and Ormerod, J. T. and Speed, T. P. and Yang, P. and Yang, J. Y. H. ",
- Title="{sc{M}erge leverages factor analysis, stable expression, and pseudoreplication to merge multiple single-cell {R}{N}{A}-seq datasets}",
- Journal="Proc. Natl. Acad. Sci. U.S.A.",
- Year="2019",
- Volume="116",
- Number="20",
- Pages="9775--9784",
- Month="May"
-}
-
-@Article{butler2018integrating,
- Author="Butler, A. and Hoffman, P. and Smibert, P. and Papalexi, E. and Satija, R. ",
- Title="{{I}ntegrating single-cell transcriptomic data across different conditions, technologies, and species}",
- Journal="Nat. Biotechnol.",
- Year="2018",
- Volume="36",
- Number="5",
- Pages="411--420",
- Month="06"
-}
-
-@ARTICLE{mcInnes2018umap,
- author = {{McInnes}, Leland and {Healy}, John and {Melville}, James},
- title = "{UMAP: Uniform Manifold Approximation and Projection for Dimension Reduction}",
- journal = {arXiv e-prints},
- keywords = {Statistics - Machine Learning, Computer Science - Computational Geometry, Computer Science - Machine Learning},
- year = "2018",
- month = "Feb",
- eid = {arXiv:1802.03426},
- pages = {arXiv:1802.03426},
-archivePrefix = {arXiv},
- eprint = {1802.03426},
- primaryClass = {stat.ML},
- adsurl = {https://ui.adsabs.harvard.edu/abs/2018arXiv180203426M},
- adsnote = {Provided by the SAO/NASA Astrophysics Data System}
-}
-
-@article{szubert2019ivis,
- author = {Szubert, B. and Cole, J. E. and Monaco, C. and Drozdov, I.},
- title = {Structure-preserving visualisation of high dimensional single-cell datasets},
- journal = {Sci Rep},
- volume = {9},
- number = {1},
- pages = {8914},
- ISSN = {2045-2322 (Electronic)
-2045-2322 (Linking)},
- DOI = {10.1038/s41598-019-45301-0},
- url = {https://www.ncbi.nlm.nih.gov/pubmed/31222035},
- year = {2019},
- type = {Journal Article}
-}
-
-@Article{scialdone2016resolving,
- Author="Scialdone, A. and Tanaka, Y. and Jawaid, W. and Moignard, V. and Wilson, N. K. and Macaulay, I. C. and Marioni, J. C. and Gottgens, B. ",
- Title="{{R}esolving early mesoderm diversification through single-cell expression profiling}",
- Journal="Nature",
- Year="2016",
- Volume="535",
- Number="7611",
- Pages="289--293",
- Month="07"
-}
-
-@article{pijuansala2019single,
- Author="Pijuan-Sala, B. and Griffiths, J. A. and Guibentif, C. and Hiscock, T. W. and Jawaid, W. and Calero-Nieto, F. J. and Mulas, C. and Ibarra-Soria, X. and Tyser, R. C. V. and Ho, D. L. L. and Reik, W. and Srinivas, S. and Simons, B. D. and Nichols, J. and Marioni, J. C. and G{\"o}ttgens, B.",
- Title="A single-cell molecular map of mouse gastrulation and early organogenesis",
- Journal="Nature",
- Year="2019",
- Volume="566",
- Number="7745",
- Pages="490--495",
- Month="02"
-}
-
-@Article{chen2016reads,
- Author="Chen, Y. and Lun, A. T. and Smyth, G. K. ",
- Title="{{F}rom reads to genes to pathways: differential expression analysis of {R}{N}{A}-{S}eq experiments using {R}subread and the edge{R} quasi-likelihood pipeline}",
- Journal="F1000Res",
- Year="2016",
- Volume="5",
- Pages="1438"
-}
-
-@Article{phipson2016robust,
- Author="Phipson, B. and Lee, S. and Majewski, I. J. and Alexander, W. S. and Smyth, G. K.",
- year = 2016,
- title = "Robust hyperparameter estimation protects against hypervariable genes and improves power to detect differential expression",
- Journal="Ann. Appl. Stat.",
- Volume=10,
- Number=2,
- Pages="946-963"
-}
-
-@Article{finak2014opencyto,
- Author="Finak, G. and Frelinger, J. and Jiang, W. and Newell, E. W. and Ramey, J. and Davis, M. M. and Kalams, S. A. and De Rosa, S. C. and Gottardo, R. ",
- Title="{{O}pen{C}yto: an open source infrastructure for scalable, robust, reproducible, and automated, end-to-end flow cytometry data analysis}",
- Journal="PLoS Comput. Biol.",
- Year="2014",
- Volume="10",
- Number="8",
- Pages="e1003806",
- Month="Aug"
-}
-
-@Article{lun2017testing,
- Author="Lun, A. T. L. and Richard, A. C. and Marioni, J. C. ",
- Title="{{T}esting for differential abundance in mass cytometry data}",
- Journal="Nat. Methods",
- Year="2017",
- Volume="14",
- Number="7",
- Pages="707--709",
- Month="Jul"
-}
-
-@article{mereu2019benchmarking,
- author = {Mereu, Elisabetta and Lafzi, Atefeh and Moutinho, Catia and Ziegenhain, Christoph and MacCarthy, Davis J. and Alvarez, Adrian and Batlle, Eduard and Sagar, and Gr{\"u}n, Dominic and Lau, Julia K. and Boutet, St{\'e}phane C. and Sanada, Chad and Ooi, Aik and Jones, Robert C. and Kaihara, Kelly and Brampton, Chris and Talaga, Yasha and Sasagawa, Yohei and Tanaka, Kaori and Hayashi, Tetsutaro and Nikaido, Itoshi and Fischer, Cornelius and Sauer, Sascha and Trefzer, Timo and Conrad, Christian and Adiconis, Xian and Nguyen, Lan T. and Regev, Aviv and Levin, Joshua Z. and Parekh, Swati and Janjic, Aleksandar and Wange, Lucas E. and Bagnoli, Johannes W. and Enard, Wolfgang and Gut, Marta and Sandberg, Rickard and Gut, Ivo and Stegle, Oliver and Heyn, Holger},
- title = {Benchmarking single-cell RNA sequencing protocols for Cell Atlas projects},
- elocation-id = {630087},
- year = {2019},
- doi = {10.1101/630087},
- journal = {bioRxiv}
-}
-
-@Article{ziegenhain2017comparative,
- Author="Ziegenhain, C. and Vieth, B. and Parekh, S. and Reinius, B. and Guillaumet-Adkins, A. and Smets, M. and Leonhardt, H. and Heyn, H. and Hellmann, I. and Enard, W. ",
- Title="{{C}omparative {A}nalysis of {S}ingle-{C}ell {R}{N}{A} {S}equencing {M}ethods}",
- Journal="Mol. Cell",
- Year="2017",
- Volume="65",
- Number="4",
- Pages="631--643",
- Month="Feb"
-}
-
-@Article{langfelder2007eigengene,
- Author="Langfelder, P. and Horvath, S. ",
- Title="{{E}igengene networks for studying the relationships between co-expression modules}",
- Journal="BMC Syst Biol",
- Year="2007",
- Volume="1",
- Pages="54",
- Month="Nov"
-}
-
-@Article{aran2019reference,
- Author="Aran, D. and Looney, A. P. and Liu, L. and Wu, E. and Fong, V. and Hsu, A. and Chak, S. and Naikawadi, R. P. and Wolters, P. J. and Abate, A. R. and Butte, A. J. and Bhattacharya, M. ",
- Title="{{R}eference-based analysis of lung single-cell sequencing reveals a transitional profibrotic macrophage}",
- Journal="Nat. Immunol.",
- Year="2019",
- Volume="20",
- Number="2",
- Pages="163--172",
- Month="02"
-}
-
-@Article{martens2013blueprint,
- Author="Martens, J. H. and Stunnenberg, H. G. ",
- Title="{{B}{L}{U}{E}{P}{R}{I}{N}{T}: mapping human blood cell epigenomes}",
- Journal="Haematologica",
- Year="2013",
- Volume="98",
- Number="10",
- Pages="1487--1489",
- Month="Oct"
-}
-
-@Article{encode2012integrated,
- Author="{The ENCODE Project Consortium}",
- Title="{{A}n integrated encyclopedia of {D}{N}{A} elements in the human genome}",
- Journal="Nature",
- Year="2012",
- Volume="489",
- Number="7414",
- Pages="57--74",
- Month="Sep"
-}
-
-@article{young2018soupx,
- author = {Young, M. D. and Behjati, S.},
- title = {Soup{X} removes ambient {R}{N}{A} contamination from droplet based single cell {R}{N}{A} sequencing data},
- year = {2018},
- journal = {bioRxiv}
-}
-
-@Article{messmer2019transcriptional,
- Author="Messmer, T. and von Meyenn, F. and Savino, A. and Santos, F. and Mohammed, H. and Lun, A. T. L. and Marioni, J. C. and Reik, W. ",
- Title="{{T}ranscriptional heterogeneity in naive and primed human pluripotent stem cells at single-cell resolution}",
- Journal="Cell Rep",
- Year="2019",
- Volume="26",
- Number="4",
- Pages="815--824",
- Month="Jan"
-}
-
-@Article{godec2016compendium,
- Author="Godec, J. and Tan, Y. and Liberzon, A. and Tamayo, P. and Bhattacharya, S. and Butte, A. J. and Mesirov, J. P. and Haining, W. N. ",
- Title="{{C}ompendium of {I}mmune {S}ignatures {I}dentifies {C}onserved and {S}pecies-{S}pecific {B}iology in {R}esponse to {I}nflammation}",
- Journal="Immunity",
- Year="2016",
- Volume="44",
- Number="1",
- Pages="194--206",
- Month="Jan"
-}
-
-@Article{tian2019benchmarking,
- Author="Tian, L. and Dong, X. and Freytag, S. and Le Cao, K. A. and Su, S. and JalalAbadi, A. and Amann-Zalcenstein, D. and Weber, T. S. and Seidi, A. and Jabbari, J. S. and Naik, S. H. and Ritchie, M. E. ",
- Title="{{B}enchmarking single cell {R}{N}{A}-sequencing analysis pipelines using mixture control experiments}",
- Journal="Nat. Methods",
- Year="2019",
- Volume="16",
- Number="6",
- Pages="479--487",
- Month="06"
-}
-
-@Article{stoeckius2017simultaneous,
- Author="Stoeckius, M. and Hafemeister, C. and Stephenson, W. and Houck-Loomis, B. and Chattopadhyay, P. K. and Swerdlow, H. and Satija, R. and Smibert, P. ",
- Title="{{S}imultaneous epitope and transcriptome measurement in single cells}",
- Journal="Nat. Methods",
- Year="2017",
- Volume="14",
- Number="9",
- Pages="865--868",
- Month="Sep"
-}
-
-@Article{lun2016csaw,
- Author="Lun, A. T. and Smyth, G. K. ",
- Title="{csaw: a {B}ioconductor package for differential binding analysis of {C}h{I}{P}-seq data using sliding windows}",
- Journal="Nucleic Acids Res.",
- Year="2016",
- Volume="44",
- Number="5",
- Pages="e45",
- Month="Mar"
-}
-
-@Article{fehlings2018multiplex,
- Author="Fehlings, M. and Chakarov, S. and Simoni, Y. and Sivasankar, B. and Ginhoux, F. and Newell, E. W. ",
- Title="{{M}ultiplex peptide-{M}{H}{C} tetramer staining using mass cytometry for deep analysis of the influenza-specific {T}-cell response in mice}",
- Journal="J. Immunol. Methods",
- Year="2018",
- Volume="453",
- Pages="30--36",
- Month="02"
-}
-
-@Article{pauken2015overcoming,
- Author="Pauken, K. E. and Wherry, E. J. ",
- Title="{{O}vercoming {T} cell exhaustion in infection and cancer}",
- Journal="Trends Immunol.",
- Year="2015",
- Volume="36",
- Number="4",
- Pages="265--276",
- Month="Apr"
-}
-
-@Article{georgiou2014promise,
- Author="Georgiou, G. and Ippolito, G. C. and Beausang, J. and Busse, C. E. and Wardemann, H. and Quake, S. R. ",
- Title="{{T}he promise and challenge of high-throughput sequencing of the antibody repertoire}",
- Journal="Nat. Biotechnol.",
- Year="2014",
- Volume="32",
- Number="2",
- Pages="158--168",
- Month="Feb"
-}
-
-@Article{rosati2017overview,
- Author="Rosati, E. and Dowds, C. M. and Liaskou, E. and Henriksen, E. K. K. and Karlsen, T. H. and Franke, A. ",
- Title="{{O}verview of methodologies for {T}-cell receptor repertoire analysis}",
- Journal="BMC Biotechnol.",
- Year="2017",
- Volume="17",
- Number="1",
- Pages="61",
- Month="07"
-}
-
-@article{luxburg2010clustering,
- title={Clustering stability: an overview},
- author={Von Luxburg, U.},
- journal={Foundations and Trends in Machine Learning},
- volume={2},
- number={3},
- pages={235--274},
- year={2010}
-}
-
-@Article{shugay2018vdjdb,
- Author="Shugay, M. and Bagaev, D. V. and Zvyagin, I. V. and Vroomans, R. M. and Crawford, J. C. and Dolton, G. and Komech, E. A. and Sycheva, A. L. and Koneva, A. E. and Egorov, E. S. and Eliseev, A. V. and Van Dyk, E. and Dash, P. and Attaf, M. and Rius, C. and Ladell, K. and McLaren, J. E. and Matthews, K. K. and Clemens, E. B. and Douek, D. C. and Luciani, F. and van Baarle, D. and Kedzierska, K. and Kesmir, C. and Thomas, P. G. and Price, D. A. and Sewell, A. K. and Chudakov, D. M. ",
- Title="{{V}{D}{J}db: a curated database of {T}-cell receptor sequences with known antigen specificity}",
- Journal="Nucleic Acids Res.",
- Year="2018",
- Volume="46",
- Number="D1",
- Pages="D419-D427",
- Month="01"
-}
-
-@Article{nestorowa2016singlecell,
- Author="Nestorowa, S. and Hamey, F. K. and Pijuan Sala, B. and Diamanti, E. and Shepherd, M. and Laurenti, E. and Wilson, N. K. and Kent, D. G. and Gottgens, B. ",
- Title="{{A} single-cell resolution map of mouse hematopoietic stem and progenitor cell differentiation}",
- Journal="Blood",
- Year="2016",
- Volume="128",
- Number="8",
- Pages="20--31",
- Month="08"
-}
-
-@Article{grun2016denovo,
- Author="Grun, D. and Muraro, M. J. and Boisset, J. C. and Wiebrands, K. and Lyubimova, A. and Dharmadhikari, G. and van den Born, M. and van Es, J. and Jansen, E. and Clevers, H. and de Koning, E. J. P. and van Oudenaarden, A. ",
- Title="{{D}e {N}ovo {P}rediction of {S}tem {C}ell {I}dentity using {S}ingle-{C}ell {T}ranscriptome {D}ata}",
- Journal="Cell Stem Cell",
- Year="2016",
- Volume="19",
- Number="2",
- Pages="266--277",
- Month="08"
-}
-
-@Article{paul2015transcriptional,
- Author="Paul, F. and Arkin, Y. and Giladi, A. and Jaitin, D. A. and Kenigsberg, E. and Keren-Shaul, H. and Winter, D. and Lara-Astiaso, D. and Gury, M. and Weiner, A. and David, E. and Cohen, N. and Lauridsen, F. K. and Haas, S. and Schlitzer, A. and Mildner, A. and Ginhoux, F. and Jung, S. and Trumpp, A. and Porse, B. T. and Tanay, A. and Amit, I. ",
- Title="{{T}ranscriptional {H}eterogeneity and {L}ineage {C}ommitment in {M}yeloid {P}rogenitors}",
- Journal="Cell",
- Year="2015",
- Volume="163",
- Number="7",
- Pages="1663--1677",
- Month="Dec"
-}
-
-@article{crowell2019discovery,
- author = {Crowell, H. L. and Soneson, C. and Germain, P.-L. and Calini, D. and Collin, L. and Raposo, C. and Malhotra, D. and Robinson, M. D.},
- title = {On the discovery of population-specific state transitions from multi-sample multi-condition single-cell RNA sequencing data},
- year = {2019},
- doi = {10.1101/713412},
- journal = {bioRxiv}
-}
-
-@Article{shao2017robust,
- Author="Shao, C. and H\"ofer, T. ",
- Title="{{R}obust classification of single-cell transcriptome data by nonnegative matrix factorization}",
- Journal="Bioinformatics",
- Year="2017",
- Volume="33",
- Number="2",
- Pages="235--242",
- Month="01"
-}
-
-@Article{kotliar2019identifying,
- Author="Kotliar, D. and Veres, A. and Nagy, M. A. and Tabrizi, S. and Hodis, E. and Melton, D. A. and Sabeti, P. C. ",
- Title="{{I}dentifying gene expression programs of cell-type identity and cellular activity with single-cell {R}{N}{A}-{S}eq}",
- Journal="Elife",
- Year="2019",
- Volume="8",
- Month="07"
-}
-
-@Article{johnstone2009consistency,
- Author="Johnstone, I. M. and Lu, A. Y. ",
- Title="{{O}n {C}onsistency and {S}parsity for {P}rincipal {C}omponents {A}nalysis in {H}igh {D}imensions}",
- Journal="J Am Stat Assoc",
- Year="2009",
- Volume="104",
- Number="486",
- Pages="682--693",
- Month="06"
-}
-
-@article{hastie1989principal,
- author="Hastie , T. and Stuetzle, W.",
- title="Principal Curves",
- journal="J Am Stat Assoc",
- volume = {84},
- number = {406},
- pages = {502-516},
- year = {1989}
-}
-
-@Article{street2018slingshot,
- Author="Street, K. and Risso, D. and Fletcher, R. B. and Das, D. and Ngai, J. and Yosef, N. and Purdom, E. and Dudoit, S. ",
- Title="{{S}lingshot: cell lineage and pseudotime inference for single-cell transcriptomics}",
- Journal="BMC Genomics",
- Year="2018",
- Volume="19",
- Number="1",
- Pages="477",
- Month="Jun"
-}
-
-@Article{sealens2019comparison,
- Author="Saelens, W. and Cannoodt, R. and Todorov, H. and Saeys, Y. ",
- Title="{{A} comparison of single-cell trajectory inference methods}",
- Journal="Nat. Biotechnol.",
- Year="2019",
- Volume="37",
- Number="5",
- Pages="547--554",
- Month="05"
-}
-
-@Article{gulati2020singlecell,
- Author="Gulati, G. S. and Sikandar, S. S. and Wesche, D. J. and Manjunath, A. and Bharadwaj, A. and Berger, M. J. and Ilagan, F. and Kuo, A. H. and Hsieh, R. W. and Cai, S. and Zabala, M. and Scheeren, F. A. and Lobo, N. A. and Qian, D. and Yu, F. B. and Dirbas, F. M. and Clarke, M. F. and Newman, A. M. ",
- Title="{{S}ingle-cell transcriptional diversity is a hallmark of developmental potential}",
- Journal="Science",
- Year="2020",
- Volume="367",
- Number="6476",
- Pages="405--411",
- Month="01"
-}
-
-@Article{guo2017slice,
- Author="Guo, M. and Bao, E. L. and Wagner, M. and Whitsett, J. A. and Xu, Y. ",
- Title="{{S}{L}{I}{C}{E}: determining cell differentiation and lineage based on single cell entropy}",
- Journal="Nucleic Acids Res.",
- Year="2017",
- Volume="45",
- Number="7",
- Pages="e54",
- Month="04"
-}
-
-@Article{teschendorff2017singlecell,
- Author="Teschendorff, A. E. and Enver, T. ",
- Title="{{S}ingle-cell entropy for accurate estimation of differentiation potency from a cell's transcriptome}",
- Journal="Nat Commun",
- Year="2017",
- Volume="8",
- Pages="15599",
- Month="06"
-}
-
-@Article{lamanno2018rna,
- Author="La Manno, G. and Soldatov, R. and Zeisel, A. and Braun, E. and Hochgerner, H. and Petukhov, V. and Lidschreiber, K. and Kastriti, M. E. and L?nnerberg, P. and Furlan, A. and Fan, J. and Borm, L. E. and Liu, Z. and van Bruggen, D. and Guo, J. and He, X. and Barker, R. and Sundstr?m, E. and Castelo-Branco, G. and Cramer, P. and Adameyko, I. and Linnarsson, S. and Kharchenko, P. V. ",
- Title="{{R}{N}{A} velocity of single cells}",
- Journal="Nature",
- Year="2018",
- Volume="560",
- Number="7719",
- Pages="494--498",
- Month="08"
-}
-
-@article{bergen2019generalizing,
- author = {Bergen, Volker and Lange, Marius and Peidli, Stefan and Wolf, F. Alexander and Theis, Fabian J.},
- title = {Generalizing RNA velocity to transient cell states through dynamical modeling},
- elocation-id = {820936},
- year = {2019},
- doi = {10.1101/820936},
- publisher = {Cold Spring Harbor Laboratory},
- URL = {https://www.biorxiv.org/content/early/2019/10/29/820936},
- eprint = {https://www.biorxiv.org/content/early/2019/10/29/820936.full.pdf},
- journal = {bioRxiv}
-}
-
-@Article{shekhar2016comprehensive,
- Author="Shekhar, K. and Lapan, S. W. and Whitney, I. E. and Tran, N. M. and Macosko, E. Z. and Kowalczyk, M. and Adiconis, X. and Levin, J. Z. and Nemesh, J. and Goldman, M. and McCarroll, S. A. and Cepko, C. L. and Regev, A. and Sanes, J. R. ",
- Title="{{C}omprehensive {C}lassification of {R}etinal {B}ipolar {N}eurons by {S}ingle-{C}ell {T}ranscriptomics}",
- Journal="Cell",
- Year="2016",
- Volume="166",
- Number="5",
- Pages="1308--1323",
- Month="Aug"
-}
-
-@article{soneson2020preprocessing,
- author = {Soneson, C. and Srivastava, A. and Patro, R. and Stadler, M. B.},
- title = {Preprocessing choices affect RNA velocity results for droplet scRNA-seq data},
- year = {2020},
- doi = {10.1101/2020.03.13.990069},
- journal = {bioRxiv}
-}
-
-@Article{steinman2002cell,
- Author="Steinman, R. A. ",
- Title="{{C}ell cycle regulators and hematopoiesis}",
- Journal="Oncogene",
- Year="2002",
- Volume="21",
- Number="21",
- Pages="3403--3413",
- Month="May"
-}
-
-@Article{kozar2004mouse,
- Author="Kozar, K. and Ciemerych, M. A. and Rebel, V. I. and Shigematsu, H. and Zagozdzon, A. and Sicinska, E. and Geng, Y. and Yu, Q. and Bhattacharya, S. and Bronson, R. T. and Akashi, K. and Sicinski, P. ",
- Title="{{M}ouse development and cell proliferation in the absence of {D}-cyclins}",
- Journal="Cell",
- Year="2004",
- Volume="118",
- Number="4",
- Pages="477--491",
- Month="Aug"
-}
-
-@Article{roccio2013predicting,
- Author="Roccio, M. and Schmitter, D. and Knobloch, M. and Okawa, Y. and Sage, D. and Lutolf, M. P. ",
- Title="{{P}redicting stem cell fate changes by differential cell cycle progression patterns}",
- Journal="Development",
- Year="2013",
- Volume="140",
- Number="2",
- Pages="459--470",
- Month="Jan"
-}
-
-@Article{soufi2016cycling,
- Author="Soufi, A. and Dalton, S. ",
- Title="{{C}ycling through developmental decisions: how cell cycle dynamics control pluripotency, differentiation and reprogramming}",
- Journal="Development",
- Year="2016",
- Volume="143",
- Number="23",
- Pages="4301--4311",
- Month="12"
-}
-
-@article{buttner2019test,
- title={A test metric for assessing single-cell RNA-seq batch correction},
- author={B{\"u}ttner, Maren and Miao, Zhichao and Wolf, F Alexander and Teichmann, Sarah A and Theis, Fabian J},
- journal={Nature methods},
- volume={16},
- number={1},
- pages={43--49},
- year={2019},
- publisher={Nature Publishing Group}
-}
-
-@book{morgan2007cell,
- title={The cell cycle: principles of control},
- author={Morgan, D. O.},
- year={2007},
- publisher={New Science Press}
-}
-
-@Article{abid2018exploring,
- Author="Abid, A. and Zhang, M. J. and Bagaria, V. K. and Zou, J. ",
- Title="{{E}xploring patterns enriched in a dataset with contrastive principal component analysis}",
- Journal="Nat Commun",
- Year="2018",
- Volume="9",
- Number="1",
- Pages="2134",
- Month="05"
-}
-
-@Article{boileau2020exploring,
- Author="Boileau, P. and Hejazi, N. S. and Dudoit, S. ",
- Title="{{E}xploring high-dimensional biological data with sparse contrastive principal component analysis}",
- Journal="Bioinformatics",
- Year="2020",
- Volume="36",
- Number="11",
- Pages="3422--3430",
- Month="Jun"
-}
-
-@Article{wu2019advantages,
- Author="Wu, H. and Kirita, Y. and Donnelly, E. L. and Humphreys, B. D. ",
- Title="{{A}dvantages of {S}ingle-{N}ucleus over {S}ingle-{C}ell {R}{N}{A} {S}equencing of {A}dult {K}idney: {R}are {C}ell {T}ypes and {N}ovel {C}ell {S}tates {R}evealed in {F}ibrosis}",
- Journal="J. Am. Soc. Nephrol.",
- Year="2019",
- Volume="30",
- Number="1",
- Pages="23--32",
- Month="01"
-}
-
-@Article{bakken2018single,
- Author="Bakken, T. E. and Hodge, R. D. and Miller, J. A. and Yao, Z. and Nguyen, T. N. and Aevermann, B. and Barkan, E. and Bertagnolli, D. and Casper, T. and Dee, N. and Garren, E. and Goldy, J. and Graybuck, L. T. and Kroll, M. and Lasken, R. S. and Lathia, K. and Parry, S. and Rimorin, C. and Scheuermann, R. H. and Schork, N. J. and Shehata, S. I. and Tieu, M. and Phillips, J. W. and Bernard, A. and Smith, K. A. and Zeng, H. and Lein, E. S. and Tasic, B. ",
- Title="{{S}ingle-nucleus and single-cell transcriptomes compared in matched cortical cell types}",
- Journal="PLoS ONE",
- Year="2018",
- Volume="13",
- Number="12",
- Pages="e0209648"
-}
-
-@Article{mcginnis2019multiseq,
- Author="McGinnis, C. S. and Patterson, D. M. and Winkler, J. and Conrad, D. N. and Hein, M. Y. and Srivastava, V. and Hu, J. L. and Murrow, L. M. and Weissman, J. S. and Werb, Z. and Chow, E. D. and Gartner, Z. J. ",
- Title="{{M}{U}{L}{T}{I}-seq: sample multiplexing for single-cell {R}{N}{A} sequencing using lipid-tagged indices}",
- Journal="Nat. Methods",
- Year="2019",
- Volume="16",
- Number="7",
- Pages="619--626",
- Month="07"
-}
-
-@Article{griffiths2018detection,
- Author="Griffiths, J. A. and Richard, A. C. and Bach, K. and Lun, A. T. L. and Marioni, J. C. ",
- Title="{{D}etection and removal of barcode swapping in single-cell {R}{N}{A}-seq data}",
- Journal="Nat Commun",
- Year="2018",
- Volume="9",
- Number="1",
- Pages="2667",
- Month="07"
-}
-
-@article{sinha2017index,
- author = {Sinha, Rahul and Stanley, Geoff and Gulati, Gunsagar S. and Ezran, Camille and Travaglini, Kyle J. and Wei, Eric and Chan, Charles K.F. and Nabhan, Ahmad N. and Su, Tianying and Morganti, Rachel M. and Conley, Stephanie D. and Chaib, Hassan and Red-Horse, Kristy and Longaker, Michael T. and Snyder, Michael P. and Krasnow, Mark A. and Weissman, Irving L.},
- title = {Index switching causes {\textquotedblleft}spreading-of-signal{\textquotedblright} among multiplexed samples in Illumina HiSeq 4000 DNA sequencing},
- elocation-id = {125724},
- year = {2017},
- doi = {10.1101/125724},
- URL = {https://www.biorxiv.org/content/early/2017/04/09/125724},
- eprint = {https://www.biorxiv.org/content/early/2017/04/09/125724.full.pdf},
- journal = {bioRxiv}
-}
-
-@Article{germain2020pipecomp,
- Author="Germain, P. L. and Sonrel, A. and Robinson, M. D. ",
- Title="{pipe{C}omp, a general framework for the evaluation of computational pipelines, reveals performant single cell {R}{N}{A}-seq preprocessing tools}",
- Journal="Genome Biol.",
- Year="2020",
- Volume="21",
- Number="1",
- Pages="227",
- Month="Sep"
-}
-
-@Article{hermann2018mammalian,
- Author="Hermann, B. P. and Cheng, K. and Singh, A. and Roa-De La Cruz, L. and Mutoji, K. N. and Chen, I. C. and Gildersleeve, H. and Lehle, J. D. and Mayo, M. and Westernstr?er, B. and Law, N. C. and Oatley, M. J. and Velte, E. K. and Niedenberger, B. A. and Fritze, D. and Silber, S. and Geyer, C. B. and Oatley, J. M. and McCarrey, J. R. ",
- Title="{{T}he {M}ammalian {S}permatogenesis {S}ingle-{C}ell {T}ranscriptome, from {S}permatogonial {S}tem {C}ells to {S}permatids}",
- Journal="Cell Rep",
- Year="2018",
- Volume="25",
- Number="6",
- Pages="1650--1667",
- Month="11"
-}
-
-@Article{kotliarov2020broad,
- Author="Kotliarov, Y. and Sparks, R. and Martins, A. J. and Mul?, M. P. and Lu, Y. and Goswami, M. and Kardava, L. and Banchereau, R. and Pascual, V. and Biancotto, A. and Chen, J. and Schwartzberg, P. L. and Bansal, N. and Liu, C. C. and Cheung, F. and Moir, S. and Tsang, J. S. ",
- Title="{{B}road immune activation underlies shared set point signatures for vaccine responsiveness in healthy individuals and disease activity in patients with lupus}",
- Journal="Nat Med",
- Year="2020",
- Volume="26",
- Number="4",
- Pages="618--629",
- Month="04"
-}
-
-@Article{zhang2020determining,
- Author="Zhang, M. J. and Ntranos, V. and Tse, D. ",
- Title="{{D}etermining sequencing depth in a single-cell {R}{N}{A}-seq experiment}",
- Journal="Nat Commun",
- Year="2020",
- Volume="11",
- Number="1",
- Pages="774",
- Month="02"
-}
-
-@article{svensson2019quantifying,
- title={Quantifying the tradeoff between sequencing depth and cell number in single-cell RNA-seq},
- author={Svensson, V. and da Veiga Beltrame, E. and Pachter, L.},
- journal={bioRxiv},
- pages={762773},
- year={2019},
- publisher={Cold Spring Harbor Laboratory}
-}
-
-@Article{srivastava2019alevin,
- Author="Srivastava, A. and Malik, L. and Smith, T. and Sudbery, I. and Patro, R. ",
- Title="{{A}levin efficiently estimates accurate gene abundances from dsc{R}{N}{A}-seq data}",
- Journal="Genome Biol",
- Year="2019",
- Volume="20",
- Number="1",
- Pages="65",
- Month="03"
-}
-
-@article{bais2020scds,
- Author = {Bais, A. S. and Kostka, D.},
- Title = {scds: computational annotation of doublets in single-cell {RNA} sequencing data},
- Volume = {36},
- Journal = {Bioinformatics},
- Year = "2020",
- Number = {4},
- Pages = {1150--1158}
-}
-
-@Article{frey2007clustering,
- Author="Frey, B. J. and Dueck, D. ",
- Title="{{C}lustering by passing messages between data points}",
- Journal="Science",
- Year="2007",
- Volume="315",
- Number="5814",
- Pages="972--976",
- Month="Feb"
-}
-
diff --git a/inst/book/style.css b/inst/book/style.css
deleted file mode 100644
index 711eafd..0000000
--- a/inst/book/style.css
+++ /dev/null
@@ -1,45 +0,0 @@
-/*
-Developed and maintained by Kevin Rue-Albrecht (@kevinrue)
-*/
-
-/*
-#0092ac blue
-#c4d931 green (on blue)
-#87b13f green (on white)
-*/
-
-p.caption {
- color: #777;
- margin-top: 10px;
-}
-p code {
- white-space: inherit;
-}
-pre {
- word-break: normal;
- word-wrap: normal;
-}
-pre code {
- white-space: inherit;
-}
-
-.book .book-summary {
- background-color: #f0f0f0;
-}
-
-/*.book .book-summary ul.summary li.part>span*/
-.book .book-summary ul.summary li span {
- color: black;
-}
-
-.book .book-body .page-wrapper .page-inner section.normal h1 {
- color: black;
-}
-
-.book .book-body .page-wrapper .page-inner section.normal h2 {
- color: #87b13f;
-}
-
-.book .book-body .page-wrapper .page-inner section.normal h3 {
- color: #1a81c2;
-}
diff --git a/inst/rebook/prefix.csv b/inst/rebook/prefix.csv
deleted file mode 100644
index 5d0fc90..0000000
--- a/inst/rebook/prefix.csv
+++ /dev/null
@@ -1 +0,0 @@
-Basic
diff --git a/inst/rebook/references.csv b/inst/rebook/references.csv
deleted file mode 100644
index 4c9b660..0000000
--- a/inst/rebook/references.csv
+++ /dev/null
@@ -1,109 +0,0 @@
-id,file,text
-cell-type-annotation,cell-type-annotation.html,7
-motivation-3,cell-type-annotation.html,7.1
-assigning-cell-labels-from-reference-data,cell-type-annotation.html,7.2
-overview-2,cell-type-annotation.html,7.2.1
-using-existing-references,cell-type-annotation.html,7.2.2
-using-custom-references,cell-type-annotation.html,7.2.3
-assigning-cell-labels-from-gene-sets,cell-type-annotation.html,7.3
-assigning-cluster-labels-from-markers,cell-type-annotation.html,7.4
-computing-gene-set-activities,cell-type-annotation.html,7.5
-fig:singler-heat-pbmc,cell-type-annotation.html,7.1
-fig:singler-cluster,cell-type-annotation.html,7.2
-fig:singler-comp-pancreas,cell-type-annotation.html,7.3
-fig:auc-dist,cell-type-annotation.html,7.4
-fig:aucell-muraro-heat,cell-type-annotation.html,7.5
-fig:violin-milk,cell-type-annotation.html,7.6
-fig:lipid-synth-violin,cell-type-annotation.html,7.7
-fig:thrsp-violin,cell-type-annotation.html,7.8
-clustering,clustering.html,5
-overview-1,clustering.html,5.1
-clustering-graph,clustering.html,5.2
-background,clustering.html,5.2.1
-implementation,clustering.html,5.2.2
-adjusting-the-parameters,clustering.html,5.2.3
-vector-quantization-with-k-means,clustering.html,5.3
-background-1,clustering.html,5.3.1
-implementation-1,clustering.html,5.3.2
-in-two-step-procedures,clustering.html,5.3.3
-hierarchical-clustering,clustering.html,5.4
-background-2,clustering.html,5.4.1
-implementation-2,clustering.html,5.4.2
-in-two-step-procedures-again,clustering.html,5.4.3
-subclustering,clustering.html,5.5
-fig:tsne-clust-graph,clustering.html,5.1
-fig:pbmc-force,clustering.html,5.2
-fig:walktrap-v-others,clustering.html,5.3
-fig:tsne-clust-kmeans,clustering.html,5.4
-fig:tsne-clust-kmeans-best,clustering.html,5.5
-fig:tsne-kmeans-graph-pbmc,clustering.html,5.6
-fig:dend-416b,clustering.html,5.7
-fig:dend-cluster,clustering.html,5.8
-fig:tsne-416b,clustering.html,5.9
-fig:tsne-kmeans-hclust-pbmc,clustering.html,5.10
-fig:dendroid-pbmc,clustering.html,5.11
-fig:tsne-kmeans-affinity-pbmc,clustering.html,5.12
-fig:ccr7-dist-memory,clustering.html,5.13
-fig:cd48-memory-expression,clustering.html,5.14
-dimensionality-reduction,dimensionality-reduction.html,4
-overview,dimensionality-reduction.html,4.1
-principal-components-analysis,dimensionality-reduction.html,4.2
-choosing-the-number-of-pcs,dimensionality-reduction.html,4.3
-visualizing-the-pcs,dimensionality-reduction.html,4.4
-non-linear-methods-for-visualization,dimensionality-reduction.html,4.5
-t-stochastic-neighbor-embedding,dimensionality-reduction.html,4.5.1
-uniform-manifold-approximation-and-projection,dimensionality-reduction.html,4.5.2
-visualization-interpretation,dimensionality-reduction.html,4.5.3
-fig:zeisel-scree,dimensionality-reduction.html,4.1
-fig:zeisel-pca,dimensionality-reduction.html,4.2
-fig:zeisel-pca-multi,dimensionality-reduction.html,4.3
-fig:tsne-brain,dimensionality-reduction.html,4.4
-fig:tsne-perplexity,dimensionality-reduction.html,4.5
-fig:umap-brain,dimensionality-reduction.html,4.6
-feature-selection,feature-selection.html,3
-motivation-1,feature-selection.html,3.1
-quantifying-per-gene-variation,feature-selection.html,3.2
-sec:spikeins,feature-selection.html,3.3
-variance-batch,feature-selection.html,3.4
-hvg-selection,feature-selection.html,3.5
-feature-selection-subsetting,feature-selection.html,3.6
-fig:trend-plot-pbmc,feature-selection.html,3.1
-fig:spike-416b,feature-selection.html,3.2
-fig:tech-pbmc,feature-selection.html,3.3
-fig:blocked-fit,feature-selection.html,3.4
-marker-detection,marker-detection.html,6
-motivation-2,marker-detection.html,6.1
-scoring-markers-by-pairwise-comparisons,marker-detection.html,6.2
-effect-sizes-for-pairwise-comparisons,marker-detection.html,6.3
-summarizing-pairwise-effects,marker-detection.html,6.4
-obtaining-the-full-effects,marker-detection.html,6.5
-using-a-log-fold-change-threshold,marker-detection.html,6.6
-marker-batch,marker-detection.html,6.7
-fig:pbmc-scored-markers-1,marker-detection.html,6.1
-fig:pbmc-scored-markers-again,marker-detection.html,6.2
-fig:pbmc-grouped-ranked-heat,marker-detection.html,6.3
-fig:pbmc-grouped-ranked-heat2,marker-detection.html,6.4
-fig:pbmc-markers-5-lfc,marker-detection.html,6.5
-fig:blocked-markers-416b,marker-detection.html,6.6
-normalization,normalization.html,2
-motivation,normalization.html,2.1
-library-size-normalization,normalization.html,2.2
-normalization-by-deconvolution,normalization.html,2.3
-spike-norm,normalization.html,2.4
-normalization-transformation,normalization.html,2.5
-fig:histlib,normalization.html,2.1
-fig:deconv-zeisel,normalization.html,2.2
-fig:norm-spike-t,normalization.html,2.3
-fig:norm-effect-malat,normalization.html,2.4
-quality-control,quality-control.html,1
-quality-control-motivation,quality-control.html,1.1
-common-choices-of-qc-metrics,quality-control.html,1.2
-identifying-low-quality-cells,quality-control.html,1.3
-fixed-qc,quality-control.html,1.3.1
-quality-control-outlier,quality-control.html,1.3.2
-other-approaches,quality-control.html,1.3.3
-quality-control-plots,quality-control.html,1.4
-quality-control-discarded,quality-control.html,1.5
-fig:qc-dist-416b,quality-control.html,1.1
-fig:qc-mito-zeisel,quality-control.html,1.2
-fig:qc-mito-spike-zeisel,quality-control.html,1.3
diff --git a/vignettes/Makefile b/vignettes/Makefile
index 16adc7c..05debb1 100644
--- a/vignettes/Makefile
+++ b/vignettes/Makefile
@@ -1,4 +1,4 @@
all: compiled
compiled:
- "${R_HOME}/bin/R" -e "work.dir <- rebook::bookCache('OSCA.basic'); handle <- rebook::preCompileBook('../inst/book', work.dir=work.dir, desc='../DESCRIPTION'); old.dir <- setwd(work.dir); bookdown::render_book('index.Rmd'); setwd(old.dir); rebook::postCompileBook(work.dir=work.dir, final.dir='../inst/doc/book', handle=handle)"
+ "${R_HOME}/bin/R" -e "final.dir <- '../inst/doc/book'; dir.create(final.dir, recursive=TRUE, showWarnings=FALSE); all.files <- list.files('../inst/book'); file.copy(file.path('../inst/book', all.files), file.path(final.dir, all.files));"
diff --git a/vignettes/stub.Rmd b/vignettes/stub.Rmd
index aabd4e9..51de4ba 100644
--- a/vignettes/stub.Rmd
+++ b/vignettes/stub.Rmd
@@ -8,7 +8,7 @@ vignette: >
```{r, echo=FALSE}
# Some compileable code, apparently needed for the engine
# to recognize this as Rmarkdown.
-link <- BiocStyle::Biocbook("OSCA.basic", label="link")
+link <- BiocStyle::Biocbook("scrapbook", label="link")
URL <- sub(".*\\((.+))", "\\1", link)
```