In debsin/dropClust: Single Cell Transcriptome Analysis

knitr::opts_chunk$set(
  collapse = TRUE,
  eval=TRUE,
  comment = "#>"
)

This latest version of dropClust has been re-implemented to include new features and performance improvements. This vignette uses a small data set from the 10X website (3K PBMC dataset here ) to demonstrate a standard pipeline. This vignette can be used as a tutorial as well.

Setting up

library(dropClust)
set.seed(0)

Loading data

dropClust loads UMI count expression data from three input files. The files follow the same structure as the datasets available from the 10X website, i.e.:

• count matrix file in sparse format
• transcriptome identifiers as a TSV file and
• gene identifiers as a TSV file

# Load Data
# path contains decompressed files 
sce <-readfiles(path = "C:/Projects/dropClust/data/pbmc3k/hg19/")

Pre-processing

dropClust performs pre-processing to remove poor quality cells and genes. dropClust is also equipped to mitigate batch-effects that may be present. The user does not need to provide any information regarding the source of the batch for individual transcriptomes. However, the batch-effect removal step is optional.

Cells are filtered based on the total UMI count in a cell specified by parameter min_count. Poor quality genes are removed based on the minimum number of cells min_count with expressions above a given threshold min_count.

# Filter poor quality cells.  A threshold th corresponds to the total count of a cell.
sce<-FilterCells(sce)
sce<-FilterGenes(sce)

Data normalization

Count normalization is then performed with the good quality genes only. Normalized expression values is computed on the raw count data in a SingleCellExperiment object, using the median normalized total count.

sce<-CountNormalize(sce)

Selecting highly variable genes

Further gene selection is carried out by ranking the genes based on its dispersion index.

sce<-RankGenes(sce, ngenes_keep = 1000)

Structure Preserving Sampling

Primary clustering is performed in a fast manner to estimate a gross structure of the data. Each of these clusters is then sampled to fine tune the clustering process.

sce<-Sampling(sce)

Gene selection based on PCA

Another gene selection is performed to reduce the number of dimensions. PCA is used to identify genes affecting major components.

# Find PCA top 200 genes. This may take some time.
sce<-RankPCAGenes(sce)

Perform clustering

Fine tuning the clustering process

By default best-fit, Louvain based clusters are returned. However, the user can tune the parameters to produce the desired number of clusters. The un-sampled transcriptomes are assigned cluster identifiers from among those identifiers produced from fine-tuning clustering. The post-hoc assignment can be controlled by setting the confidence value conf. High conf values will assign cluster identifiers to only those transcriptomes sharing a majority of common nearest neighbours.

sce<-Cluster(sce, method = "default", conf = 0.5)

Visualizing clusters

Compute 2D embeddings for samples followed by post-hoc clustering.

sce<-PlotEmbedding(sce, embedding = "umap", spread = 10, min_dist = 0.1)
ScatterPlot(sce,title = "Clusters")

Find cluster specific Differentially Expressed genes

DE_genes_all = FindMarkers(sce, selected_clusters=NA, lfc_th = 1, q_th =0.001, nDE=30)

write.csv(DE_genes_all$genes, 
          file = file.path(tempdir(),"ct_genes.csv"),
          quote = FALSE)

Plot hand picked marker genes

marker_genes = c("S100A8", "GNLY", "PF4")

p<-PlotMarkers(sce, marker_genes)

Heat map of top DE genes from each cluster

# Draw heatmap
#############################
p<-PlotHeatmap(sce, DE_res = DE_genes_all$DE_res,nDE = 10)

print(p)

debsin/dropClust documentation built on Nov. 4, 2019, 10:22 a.m.