Detection of genes with differential expression dispersion unravels the role of autophagy in cancer progression

Abstract

The majority of gene expression studies focus on the search for genes whose mean expression is different between two or more populations of samples in the so-called "differential expression analysis" approach. However, a difference in variance in gene expression may also be biologically and physiologically relevant. In the classical statistical model used to analyze RNA-sequencing (RNA-seq) data, the dispersion, which defines the variance, is only considered as a parameter to be estimated prior to identifying a difference in mean expression between conditions of interest. Here, we propose to evaluate four recently published methods, which detect differences in both the mean and dispersion in RNA-seq data. We thoroughly investigated the performance of these methods on simulated datasets and characterized parameter settings to reliably detect genes with a differential expression dispersion. We applied these methods to The Cancer Genome Atlas datasets. Interestingly, among the genes with an increased expression dispersion in tumors and without a change in mean expression, we identified some key cellular functions, most of which were related to catabolism and were overrepresented in most of the analyzed cancers. In particular, our results highlight autophagy, whose role in cancerogenesis is context-dependent, illustrating the potential of the differential dispersion approach to gain new insights into biological processes and to discover new biomarkers.

Fig 1. Ability to identify differentially dispersed genes.

The performances of Levene’s test, MDSeq, DiPhiSeq, GAMLSS and DiffDist for differential dispersion detection in gene expression data, as measured by the area under the ROC curve (AUC), were assessed using 10 replicates of simulated datasets composed of highly and lowly differentially expressed genes between two sample populations of equal size.

Fig 2. Ability to detect differential dispersion for lowly differentially expressed genes.

False discovery rate (FDR) and true positive rate (TPR) of Levene’s test, MDSeq, DiPhiSeq, GAMLSS and DiffDist for differential dispersion detection in simulated datasets composed of lowly differentially expressed genes between two sample populations of equal size. The performances were assessed using 10 replicates of simulated datasets per parameter setting.

Fig 3. Differentially dispersed genes correctly identified by the evaluated methods among lowly differentially expressed genes.

(A) Intersections of sets of differentially dispersed (DD) genes correctly identified by Levene’s test, MDSeq, DiPhiSeq, GAMLSS and DiffDist. (B): Correctness of dispersion log₂-fold change sign of DD genes correctly identified by the different methods. (C) Real mean and dispersion log₂-fold changes and estimated dispersion log₂-fold changes of DD genes correctly identified by GAMLSS and DiffDist. (D) Correctness of dispersion log₂-fold change signs according to Levene’s test, MDSeq and DiPhiSeq for DD genes correctly identified by GAMLSS and DiffDist with incorrect dispersion log₂-fold change sign. Simulated datasets are composed of lowly differentially expressed genes with a mean fold change of expression between 1 and 1.5 between two populations of 50 samples. Values indicated at the middle of the bars are percentages of the corresponding categories of genes over the entire sets of analyzed genes. All results relate to 10 replicates of simulated datasets, e.g. the counts and percentages are averaged over all the replicates.

Fig 4. Differentially dispersed genes among non-differentially expressed genes for each TCGA dataset.

(A) Number of differentially expressed (DE) genes separated between those upregulated in tumors (DE+) and those downregulated in tumors (DE-) detected by MDSeq per TCGA dataset. (B) Number of differentially dispersed (DD) genes among non-DE genes separated between those overdispersed in tumors (DD+) and those underdispersed in tumors (DD-), as detected by Levene’s test, MDSeq, DiPhiSeq, GAMLSS and DiffDist, per TCGA dataset.

Fig 5. Overdispersed genes in tumors identified by the evaluated methods among non-differentially expressed genes.

Intersections of sets of overdispersed genes in tumors identified by Levene’s test, MDSeq, DiPhiSeq, GAMLSS and DiffDist among non-differentially expressed genes for (A) the kidney renal clear cell carcinoma dataset (TCGA-KIRC) and (B) the kidney renal papillary cell carcinoma dataset (TCGA-KIRP). Non-differentially expressed genes were identified by MDSeq.

Fig 6. Enriched GO terms among overdispersed genes in tumors identified by the evaluated methods.

Top 40 representative enriched Gene Ontology (GO) terms among overdispersed genes in tumors (DD+) among non-differentially expressed (non-DE) genes, ordered first by the number of datasets for which they are enriched (decreasing order) and second by the mean p-values of enrichment across all datasets (increasing order). Non-DE genes were identified using MDSeq, and DD+ genes were identified among non-differentially expressed genes by at least one of the evaluated methods, i.e. Levene’s test, MDSeq, DiPhiSeq, GAMLSS and DiffDist.