SplicingFactory-splicing diversity analysis for transcriptome data

Abstract

Motivation: Alternative splicing contributes to the diversity of RNA found in biological samples. Current tools investigating patterns of alternative splicing check for coordinated changes in the expression or relative ratio of RNA isoforms where specific isoforms are up- or down-regulated in a condition. However, the molecular process of splicing is stochastic and changes in RNA isoform diversity for a gene might arise between samples or conditions. A specific condition can be dominated by a single isoform, while multiple isoforms with similar expression levels can be present in a different condition. These changes might be the result of mutations, drug treatments or differences in the cellular or tissue environment. Here, we present a tool for the characterization and analysis of RNA isoform diversity using isoform level expression measurements.

Results: We developed an R package called SplicingFactory, to calculate various RNA isoform diversity metrics, and compare them across conditions. Using the package, we tested the effect of RNA-seq quantification tools, quantification uncertainty, gene expression levels and isoform numbers on the isoform diversity calculation. We analyzed a set of CD34+ hematopoietic stem cells and myelodysplastic syndrome samples and found a set of genes whose isoform diversity change is associated with SF3B1 mutations.

Availability and implementation: The SplicingFactory package is freely available under the GPL-3.0 license from Bioconductor for the Windows, MacOS and Linux operating systems (https://www.bioconductor.org/packages/release/bioc/html/SplicingFactory.html).

Supplementary information: Supplementary data are available at Bioinformatics online.

Fig. 1.

Correlation of various diversity metrics with gene expression and isoform number. Each dot represents a gene with at least two isoforms and at least one protein coding isoform. The x axis shows the log10(mean TPM) expression across the 17 control samples from the SRP133442 dataset, while the y axis shows a given mean diversity metric across the same samples. The color of the dots shows the log₂(isoform number) for the gene. Gene with >20 isoforms were assigned the value 20. The purple line shows the smoothed conditional mean of the data using a generalized additive model. Naive—nonnormalized naive-entropy, Laplace—nonnormalized Laplace-entropy, Gini—Gini-index, Simpson—Simpson-index, Inv. Simpson—inverse Simpson-index, Naive norm—normalized naive-entropy, Laplace norm: normalized Laplace-entropy

Fig. 2.

Various factors influencing diversity metrics. (A) Consistency of diversity metrics using the Kallisto, Salmon and Salmon-SAF transcript quantification methods. Each boxplot shows the Spearman correlation values for a specific pairwise comparison in the 17 control samples from the SRP133442 dataset, showing two different diversity metrics. K—Kallisto, S—Salmon, SS—Salmon-SAF. (B) Influence of expression estimation uncertainty on diversity calculation. The x axis shows the expression MAD, while the y axis shows the naive normalized entropy MAD of selected marker genes in a single sample. Dot color shows the original nonbootstrap expression estimates for the gene. The corr. value on the panels is the Spearman-correlation of the expression MAD and the diversity MAD. (C) Influence of filtering low-expression transcripts before diversity calculation. The x axis shows the % criteria used for filtering transcripts (see Section 2), while the y axis shows the Spearman correlation of the diversity metrics using filtered results and the original unfiltered results. Vertical lines at the dots show bootstrap replicate confidence interval

Fig. 3.

Performance benchmarks and comparison of SplicingFactory to other tools. (A) Memory usage of SplicingFactory with increasing sample number for all diversity metrics while calculating the diversity values or also calculating differential diversity between sample groups. The x axis shows the increasing number of samples used, while the y axis shows the maximum amount of memory used. (B) Total elapsed time for calculating the diversity values or also calculating differential diversity between sample groups for all diversity metrics. The x axis shows the increasing number of samples used, while the y axis shows the total elapsed time for the calculation. (C) Spearman correlation of SpliceHetero entropy values with three different diversity metrics and Salmon expression estimates for the 17 MDS samples from the SRP133442 dataset. (D) Spearman-correlation of average modified variance values (E1 and E2) as calculated by the GSReg.SEVA function and average diversity values for three different diversity metrics for the 17 control and MDS samples from the SRP133442 dataset. (E) Spearman correlation of Whippet entropy and average diversity values for three different diversity metrics for the 17 control and MDS samples from the SRP133442 dataset

Fig. 4.

Differential diversity analysis results using MDS data. (A) Differential diversity results comparing SF3B1 mutated and nonmutated MDS samples for three different diversity metrics and Salmon expression estimates. The x axis shows the mean diversity across all samples, while the y axis shows the difference of means between the SF3B1 mutated and nonmutated MDS sample groups. Each gene is a single gray dot, and significant changes are highlighted in red. Significant changes are defined as |mean difference| > 0.1 and Wilcoxon-test adjusted P-value <0.05. (B) Enrichment of the Hay et al. (2018) marker gene sets in the differential diversity results separated by increasing or decreasing diversity between the SF3B1 mutated and nonmutated sample groups. The x axis shows the different marker gene sets, while the y axis shows the % of significantly changing genes that are falling into a specific marker gene set. Significant enrichment of a set is marked with an asterisk (Benjamini-Hochberg adjusted Fisher-test P-value <0.05)