DoubletFinder

repository·master·Indexed 20 days ago

https://github.com/chris-mcginnis-ucsf/doubletfinder

An R package designed to predict and remove doublet cells from single-cell RNA sequencing data. It integrates with the Seurat ecosystem (compatible with Seurat v5 and Seurat >= 2.0) and utilizes a pN-pK parameter sweep strategy to identify optimal parameters for doublet detection. The tool employs a four-step process involving the generation of artificial doublets, pre-processing, pANN computation via PCA, and thresholding to classify cells as singlets or doublets.

Tokens
3.2K
Snippets
5
Records
11
Agent score
19%

What's inside DoubletFinder

  1. How DoubletFinder works

    master

    DoubletFinder predicts doublets in single-cell RNA sequencing data through a four-step process:

    1. Generate artificial doublets: Creates synthetic doublets from the existing scRNA-seq data.
    2. Pre-process: Merges the real and artificial data for analysis.
    3. Compute pANN: Performs PCA and uses the PC distance matrix to calculate each cell's proportion of artificial k nearest neighbors (pANN).
    4. Thresholding: Ranks the pANN values and applies a threshold based on the expected number of doublets to make final predictions.
  2. How to select the optimal pK value

    master

    The pK parameter (PC neighborhood size) is critical for accuracy. Because real-world datasets lack ground-truth doublet labels, DoubletFinder uses the mean-variance-normalized bimodality coefficient (BCmvn) to find the optimal pK.

    The Strategy:

    1. Perform a pN-pK parameter sweep.
    2. Calculate the BCmvn for the results.
    3. Select the pK value that corresponds to the single, easily-discernible maximum in the BCmvn distribution. This maximum typically coincides with the pK values that maximize AUC in validated datasets.
  3. Best practices for input scRNA-seq data

    master

    To ensure accurate doublet prediction, follow these data preparation guidelines:

    Data Aggregation

    • Do not apply DoubletFinder to aggregated data representing multiple distinct samples (e.g., different 10X lanes from different cell lines). This creates artificial doublets that cannot exist in reality and skews results.
    • It is acceptable to run DoubletFinder on data from a single sample that was split across multiple 10X lanes.
    • Avoid running DoubletFinder on integrated Seurat objects.

    Quality Control

    Ensure input data is cleared of low-quality cell clusters. A recommended workflow is:

    1. Manually threshold raw gene expression matrices based on RNA nUMIs.
    2. Pre-process using standard workflows.
    3. Identify and remove clusters characterized by low RNA UMIs, high % mitochondrial reads, or uninformative marker genes.
    4. Re-process the cleaned data before running DoubletFinder.
  4. Identify optimal pK using paramSweep, summarizeSweep, and find.pK

    master

    To identify the optimal pK value for doublet detection, you must perform a parameter sweep. This process can be done either without ground-truth data or by using existing ground-truth labels (e.g., from sample multiplexing or in silico geneotyping).

    Without Ground-Truth

    Use paramSweep to generate results, summarizeSweep to aggregate them, and find.pK to identify the best parameter.

    With Ground-Truth

    If you have ground-truth labels (a vector containing "Singlet" and "Doublet" calls), pass them to summarizeSweep using the GT = TRUE and GT.calls arguments to improve the accuracy of the pK identification.

    ## pK Identification (no ground-truth) ---------------------------------------------------------------------------------------
    sweep.res.list_kidney <- paramSweep(seu_kidney, PCs = 1:10, sct = FALSE)
    sweep.stats_kidney <- summarizeSweep(sweep.res.list_kidney, GT = FALSE)
    bcmvn_kidney <- find.pK(sweep.stats_kidney)
    
    ## pK Identification (ground-truth) ------------------------------------------------------------------------------------------
    sweep.res.list_kidney <- paramSweep(seu_kidney, PCs = 1:10, sct = FALSE)
    gt.calls <- seu_kidney@meta.data[rownames(sweep.res.list_kidney[[1]]), "GT"]
    sweep.stats_kidney <- summarizeSweep(sweep.res.list_kidney, GT = TRUE, GT.calls = gt.calls)
    bcmvn_kidney <- find.pK(sweep.stats_kidney)
  5. Estimating the expected number of doublets (nExp)

    master

    DoubletFinder is sensitive to heterotypic doublets (different cell states) but insensitive to homotypic doublets (similar cell states).

    Because Poisson statistical estimates (based on cell loading density) are agnostic to homotypic doublets, they tend to overestimate the number of detectable doublets.

    Recommendation: To get a more realistic estimate, use literature-supported cell type annotations to model the proportion of homotypic doublets in your data. Using Poisson estimates with and without homotypic adjustment provides a 'bookend' range for the real detectable doublet rate.

  6. Complete DoubletFinder workflow example

    master

    This example demonstrates a full workflow: pre-processing a Seurat object (using standard or SCTransform methods), performing a pK parameter sweep, estimating homotypic proportions, and running doubletFinder with adjusted doublet counts.

    ## Pre-process Seurat object (standard) --------------------------------------------------------------------------------------
    seu_kidney <- CreateSeuratObject(kidney.data)
    seu_kidney <- NormalizeData(seu_kidney)
    seu_kidney <- FindVariableFeatures(seu_kidney, selection.method = "vst", nfeatures = 2000)
    seu_kidney <- ScaleData(seu_kidney)
    seu_kidney <- RunPCA(seu_kidney)
    seu_kidney <- RunUMAP(seu_kidney, dims = 1:10)
    
    ## Pre-process Seurat object (sctransform) -----------------------------------------------------------------------------------
    seu_kidney <- CreateSeuratObject(kidney.data)
    seu_kidney <- SCTransform(seu_kidney)
    seu_kidney <- RunPCA(seu_kidney)
    seu_kidney <- RunUMAP(seu_kidney, dims = 1:10)
    
    ## pK Identification (no ground-truth) ---------------------------------------------------------------------------------------
    sweep.res.list_kidney <- paramSweep(seu_kidney, PCs = 1:10, sct = FALSE)
    sweep.stats_kidney <- summarizeSweep(sweep.res.list_kidney, GT = FALSE)
    bcmvn_kidney <- find.pK(sweep.stats_kidney)
    
    ## pK Identification (ground-truth) ------------------------------------------------------------------------------------------
    sweep.res.list_kidney <- paramSweep(seu_kidney, PCs = 1:10, sct = FALSE)
    gt.calls <- seu_kidney@meta.data[rownames(sweep.res.list_kidney[[1]]), "GT"]
    sweep.stats_kidney <- summarizeSweep(sweep.res.list_kidney, GT = TRUE, GT.calls = gt.calls)
    bcmvn_kidney <- find.pK(sweep.stats_kidney)
    
    ## Homotypic Doublet Proportion Estimate -------------------------------------------------------------------------------------
    homotypic.prop <- modelHomotypic(annotations)
    nExp_poi <- round(0.075*nrow(seu_kidney@meta.data))
    nExp_poi.adj <- round(nExp_poi*(1-homotypic.prop))
    
    ## Run DoubletFinder with varying classification stringencies ----------------------------------------------------------------
    seu_kidney <- doubletFinder(seu_kidney, PCs = 1:10, pN = 0.25, pK = 0.09, nExp = nExp_poi, reuse.pANN = NULL, sct = FALSE)
    seu_kidney <- doubletFinder(seu_kidney, PCs = 1:10, pN = 0.25, pK = 0.09, nExp = nExp_poi.adj, reuse.pANN = "pANN_0.25_0.09_913", sct = FALSE)
  7. Run DoubletFinder on a Seurat object

    master

    The doubletFinder function performs the final doublet classification. It requires a Seurat object and several parameters:

    • PCs: The principal components to use (e.g., 1:10).
    • pN: The probability of a doublet being a singlet.
    • pK: The optimal pK value identified during the parameter sweep.
    • nExp: The expected number of doublets. This can be the raw estimate or an adjusted estimate (nExp_poi.adj) that accounts for homotypic doublets.
    • reuse.pANN: Allows reusing a previously calculated pANN (probability of being an Artificial Neural Network doublet) string to maintain consistency across different classification stringencies.
    • sct: Boolean indicating if SCTransform normalization was used.
    ## Run DoubletFinder with varying classification stringencies ----------------------------------------------------------------
    seu_kidney <- doubletFinder(seu_kidney, PCs = 1:10, pN = 0.25, pK = 0.09, nExp = nExp_poi, reuse.pANN = NULL, sct = FALSE)
    seu_kidney <- doubletFinder(seu_kidney, PCs = 1:10, pN = 0.25, pK = 0.09, nExp = nExp_poi.adj, reuse.pANN = "pANN_0.25_0.09_913", sct = FALSE)
  8. Estimate homotypic doublet proportion with modelHomotypic

    master

    The modelHomotypic function estimates the proportion of homotypic doublets (doublets consisting of two cells of the same type) based on cell type annotations. This estimate is used to adjust the expected number of doublets to avoid over-correcting for cells that are naturally similar.

    ## Homotypic Doublet Proportion Estimate -------------------------------------------------------------------------------------
    homotypic.prop <- modelHomotypic(annotations)           ## ex: annotations <- seu_kidney@meta.data$ClusteringResults
  9. DoubletFinder input arguments

    master

    When using DoubletFinder, you will need to provide or configure the following arguments:

    • seu: A fully-processed Seurat object. It must have already undergone NormalizeData, FindVariableGenes, ScaleData, RunPCA, and RunTSNE.
    • PCs: The number of statistically-significant principal components, specified as a range (e.g., PCs = 1:10).
    • pN: The proportion of generated artificial doublets relative to the merged real-artificial data. The default is 25%.
    • pK: The PC neighborhood size used to compute pANN, expressed as a proportion of the merged real-artificial data. This value should be adjusted for each dataset using the pN-pK parameter sweep strategy.
    • nExp: The pANN threshold used for final doublet/singlet predictions. This is typically estimated from cell loading densities and adjusted for homotypic doublets.