scCoExpress: A Sparsity-Aware R Package for Gene Co-Expression Analysis in Single-Cell Datasets

Implementation

Broadly, the scCoExpress algorithm works by modeling a series of pairwise null interaction strength distributions through permutations, and scaling the interaction strength of target gene pairs against a depth-appropriate background distribution. For example, the interaction strength of a pair of genes which are both expressed in roughly 10% of cells would be scaled against a null distribution of randomly selected gene pairs with abundances near 10%. This establishes a background rate at which random genes at a given pair of abundances co-express, and can be used to assess the degree to which a given gene pair co-localizes relative to chance at a particular sparsity level. Co-expression intensity for all gene comparisons, both the null and target comparisons, are quantified with a modified version of the Manders’ Overlap Coefficient⁵ originally developed for assessing fluorophore co-localization in confocal arrays:

Where c represents all cells in the experiment, A represents the vector of normalized expression of gene A in all cells k, and B represents the vector of normalized expression of gene B in all cells k. For each target gene comparison, a MOC Ratio is calculated as:

A Z-score is subsequently calculated for each comparison’s MOC against the null MOC distribution. In cases where many gene pairs are queried, the Z-score is adjusted to a predicted value (“ZAdj”) based on the linear regression of all Z-scores against the raw MOC ratio, in order to reduce noise. Finally, p-values are quantified from Zadj using a standard normal table.

Background MOC distributions are tabulated in one of two ways depending on the user’s settings. The default and recommended way is via “partitions,” wherein all genes in the input are grouped into N equally-sized groups based on the fraction of cells in which each gene is detected, excluding genes with expression extremely close to or equal to either 0% or 100%. One background distribution is then calculated for each pairwise comparison of partitions, and target gene pairs are subsequently normalized against the most appropriate null set. Users also have an option to run scCoExpress in “local” mode, wherein an independent null distribution is generated for each target gene set based on the genes with the closest expression to both targets. Local mode offers improved resolution in many cases, but is substantially less memory efficient when querying a large number of gene pairs

Operation

scCoExpress is currently available for download on Github. It is a simple and user-friendly R package, which in its most minimal form can be run with a single command, and then visualized with a single subsequent command. The only required inputs are a Seurat object with a normalized assay slot, and a vector of genes:

# Install and load libraries
remotes::install_github("katlande/scCoExpress")
library(Seurat)
library(ggplot2)
library(scCoExpress)
# Run co-expression with default settings, returns a data.frame:
co <- CoExpress(obj = pbmc, target_genes = c("gene1", "gene2", "genen"))
# Visualize co-expression as a correlation plot, returns a grob:
plotCoExpr(co)

scCoExpress will query co-expression across all cells in the supplied Seurat object. To restrict the analysis to a specific cell type or variable(s), one must subset the input object as desired prior to running CoExpress(). A full vignette going over input parameters and QC has been provided on GitHub.

Limitations

Tabula Muris is composed of whole-transcriptome single-cell sequencing data from 20 mouse tissues and organs, se-quenced with either the 10x GemCode Single-Cell 3’ Gel Bead and Library V2 Kit or the SMART-seq2 protocol. Each Tabula Muris set contains between 408 and 11,269 cells. Applications of scCoExpress to other species, sequencing strategies, or dataset sizes have not been directly assessed in this article

Additionally, while it is possible to use scCoExpress to query all possible gene pairs in the genome, this is not recommended due to the memory requirements needed to bootstrap local distributions for n(n-1)/2 comparisons where n equals the total transcriptome size. The intended use case for scCoExpress is to assess pairwise co-expression relationships between a specific set of genes of interest, or to query the whole transcriptome against a specific gene of interest. For pairwise queries of extremely large gene sets, it is recommended to run partition mode with a larger number of partitions and extreme.skipping set to TRUE.