Motivation

Differential expression (DE) analysis in single-cell transcriptomics reveals cell-type-specific responses. Recent studies have raised concerns about current methods used to identify differentially expressed genes in this context.

• scRNA sequencing provides the absolute abundances of RNA molecules in single cells, but normalization - a pre-processing step inherited from the bulk RNA-seq era - reduces this information and returns data as relative abundances, which may mask important differences among cell types

• The majority of single-cell DE analysis methods are susceptible to generating false discoveries. This is mainly due to the lack of accounting for variations between biological replicates, commonly referred to as “donor effects”.

  • Batch effects are often estimated from leading principal components, representing a consensus from most genes.
  • Pseudo-bulk analysis ignores within-sample heterogeneity by treating donor effects as a fixed effect and assuming that each cell from the same donor is equally affected.

• Clustering and DE analysis are different problems. The current commonly workflow works well in clustering, but cannot guarantee success in downstream analysis.

• Excessive zeros are usually considered as “drop-outs”, while they are actually informative in cell-type heterogeneity. Ignoring zeros in single-cell gene expression data discards valuable information for any analysis.

We provide a generalized linear mixed model framework (GLMM) to detect differentially expressed genes (DEGs) between two given cell types. The model takes donor-specific variations as random effects, and uses raw UMI counts to prevent biases in DE analysis.

Past versions of panel3.png

GLMM for DE analysis

New criteria on determining DEGs

We proposed new criteria that based on the convention and also the gene mean and the difference in mean. If the log2 gene mean in two groups is lower than a certain value (-2.25 in case study 1) and the log2 mean difference is smaller than a threshold (-1 in case study 1), the genes would not be considered as a DEGs. These can also be used as a filter before any DE analysis to speed up the computation. Both of these criteria are adjustable, depending on the dataset’s performance and characteristics. An examination in heatmaps and mean difference against mean plot in advanced can be helpful to determine the thresholds when analyzing a new dataset. More details can be found here.

Benchmarked methods

In this project, we compare a few methods performing the DE analysis results. Our comparison encompassed Poisson-glmm and Binomial-glmm from the new paradigm, as well as pseudo-bulk approaches including DESeq2 and edgeR. Additionally, we assessed the performance of single cell specific tools including MAST, Wilcox in Seurat, and linear mixed models in Muscat. More details can be found here.

Analysis on case study 1

• About data

In case study 1, a 10X scRNA-seq dataset of post-menopausal fallopian tubes, with 57,182 cells sourced from five donors, covering 29,382 genes is analyzed. The 20 clusters are obtained via HIPPO algorithm. There is no pre-filtering procedure applied on this dataset, except for built-in filtering steps in each method. We use sctransform to get the VST data, and the integration workflow provided by Seurat to obtain the integrated data.

• Analyses

  • Data summary
    
  • group12_13
    
  • group2_19
    
  • group8_17&2_19
    
  • group12_19
    
  • group8_17

Analysis on case study 2

• About data In case study 2, the dataset contains 10X droplet-based scRNA-seq PBCM data from 8 Lupus patients obtained before and after 6h-treatment with IFN-$\beta$. After removing undetected genes and lowly expressed genes (less than 10 cells expressing more than 1), the dataset consists of 29065 cells and 7661 genes. The integrated data is replaced by log2-transformed normalized expression values obtained via computeLibrarayFactors and logNormCounts functions in Muscat.

• Analyses

  • Data summary
  • CD14+ Monocytes
  • B cells
  • P-values under null setting

Session information

sessionInfo()

R version 4.2.2 (2022-10-31)
Platform: x86_64-apple-darwin17.0 (64-bit)
Running under: macOS Big Sur ... 10.16

Matrix products: default
BLAS:   /Library/Frameworks/R.framework/Versions/4.2/Resources/lib/libRblas.0.dylib
LAPACK: /Library/Frameworks/R.framework/Versions/4.2/Resources/lib/libRlapack.dylib

locale:
[1] en_US.UTF-8/en_US.UTF-8/en_US.UTF-8/C/en_US.UTF-8/en_US.UTF-8

attached base packages:
[1] stats     graphics  grDevices utils     datasets  methods   base     

loaded via a namespace (and not attached):
 [1] Rcpp_1.0.11       rstudioapi_0.15.0 whisker_0.4.1     knitr_1.27       
 [5] magrittr_2.0.3    workflowr_1.7.0   R6_2.5.1          rlang_1.1.1      
 [9] fastmap_1.1.1     fansi_1.0.4       stringr_1.5.0     tools_4.2.2      
[13] xfun_0.39         png_0.1-8         utf8_1.2.3        cli_3.6.1        
[17] jquerylib_0.1.4   git2r_0.32.0      htmltools_0.5.5   rprojroot_2.0.3  
[21] yaml_2.3.7        digest_0.6.33     tibble_3.2.1      lifecycle_1.0.3  
[25] later_1.3.1       vctrs_0.6.3       sass_0.4.7        promises_1.2.0.1 
[29] fs_1.6.3          glue_1.6.2        cachem_1.0.8      evaluate_0.21    
[33] rmarkdown_2.23    stringi_1.7.12    pillar_1.9.0      compiler_4.2.2   
[37] bslib_0.5.0       jsonlite_1.8.7    httpuv_1.6.11     pkgconfig_2.0.3