Variance component score test for time-course gene set analysis of longitudinal RNA-seq data

TL;DR

This paper introduces tcgsaseq, a novel, model-free statistical method for detecting longitudinal gene set changes in RNA-seq data, effectively handling heteroscedasticity and covariates without assuming specific distributions.

Contribution

The paper presents tcgsaseq, a new variance component score test for longitudinal RNA-seq analysis that is efficient, nonparametric, and improves stability and power over existing methods.

Findings

tcgsaseq controls type I error effectively.

It shows higher power than ROAST, edgeR, and DESeq2.

The method is computationally efficient and available in R.

Abstract

As gene expression measurement technology is shifting from microarrays to sequencing, the statistical tools available for their analysis must be adapted since RNA-seq data are measured as counts. Recently, it has been proposed to tackle the count nature of these data by modeling log-count reads per million as continuous variables, using nonparametric regression to account for their inherent heteroscedasticity. Adopting such a framework, we propose tcgsaseq, a principled, model-free and efficient top-down method for detecting longitudinal changes in RNA-seq gene sets. Considering gene sets defined a priori, tcgsaseq identifies those whose expression vary over time, based on an original variance component score test accounting for both covariates and heteroscedasticity without assuming any specific parametric distribution for the transformed counts. We demonstrate that despite the presence of a nonparametric component, our test statistic has a simple form and limiting distribution, and both may be computed quickly. A permutation version of the test is additionally proposed for very small sample sizes. Applied to both simulated data and two real datasets, the proposed method is shown to exhibit very good statistical properties, with an increase in stability and power when compared to state of the art methods ROAST, edgeR and DESeq2, which can fail to control the type I error under certain realistic settings. We have made the method available for the community in the R package tcgsaseq.