SpliceTools, a suite of downstream RNA splicing analysis tools to investigate mechanisms and impact of alternative splicing

Abstract

As a fundamental aspect of normal cell signaling and disease states, there is great interest in determining alternative splicing (AS) changes in physiologic, pathologic, and pharmacologic settings. High throughput RNA sequencing and specialized software to detect AS has greatly enhanced our ability to determine transcriptome-wide splicing changes. Despite the richness of this data, deriving meaning from sometimes thousands of AS events is a substantial bottleneck for most investigators. We present SpliceTools, a suite of data processing modules that arms investigators with the ability to quickly produce summary statistics, mechanistic insights, and functional significance of AS changes through command line or through an online user interface. Utilizing RNA-seq datasets for 186 RNA binding protein knockdowns, nonsense mediated RNA decay inhibition, and pharmacologic splicing inhibition, we illustrate the utility of SpliceTools to distinguish splicing disruption from regulated transcript isoform changes, we show the broad transcriptome footprint of the pharmacologic splicing inhibitor, indisulam, we illustrate the utility in uncovering mechanistic underpinnings of splicing inhibition, we identify predicted neo-epitopes in pharmacologic splicing inhibition, and we show the impact of splicing alterations induced by indisulam on cell cycle progression. Together, SpliceTools puts rapid and easy downstream analysis at the fingertips of any investigator studying AS.

Figure 1.

SpliceTools schematic. (A) The six core skipped exon analysis tools are shown with their respective outputs indicated to the right. Tools can either be used individually or all six can be called using one command with SEMedley. (B) Listing of core tools for retained intron analysis and their outputs are indicated. RIMedley can be used to run all three retained intron tools using a single command. (C) SpliceCompare analyzes any combination of A5SS, A3SS, MXE, RI, or SE files for comparative analyses of AS profiles across datasets. A parallel processing version of SpliceCompare is provided that utilizes 5 processors to simultaneously analyze A5SS, A3SS, MXE, RI and SE files if the user has installed the Parallel::Forkmanager module.

Figure 2.

Transcriptome footprint of SE and RI for RBP knockdowns and pharmacologic splicing inhibitors. SEFractionExpressed and RIFractionExpressed were used with input parameters of minimum TPMs = 3 and an AS threshold FDR < 0.0005. Colors and intensities indicate the proportion of significantly changed AS events with increased exon skipping or intron retention. (A) fraction of expressed genes with increased exon skipping (–IncDiff (Inclusion Difference)). (B) fraction of expressed genes with increased intron retention (+IncDiff) events.

Figure 3.

FractionUnannotated. Plots indicate the fraction of statistically significant (FDR < 0.0005) increased (blue) or decreased (green) exon skipping events that are not denoted in the input annotation file.

Figure 4.

NumberSkipped. Distribution of multiple exon skipping across datasets. SENumberSkipped was run using an FDR cutoff of 0.0005. The plot includes experiments with >100 detected -IncDiff events (increase in skipping). The splicing inhibitor, indisulam, shows higher numbers of exons skipped. Increases in skipped isoform levels upon NMD inhibition (e.g. CC115 and SMG6/7 knockdown), more frequently involve single exon skipping.

Figure 5.

IntronExonSizes. Distributions of flanking intron and skipped exon sizes across ENCODE RBP knockdown dataset. SEIntronExonSizes was run using an FDR cutoff of 0.0005. Average sizes for experiments with greater than 100 -IncDiff (increase in skipping) events were plotted. (A) Correlation between upstream and downstream flanking intron sizes. Average flanking intron (B) and skipped exon (C) sizes for events with increased skipping are plotted for each knockdown. Average sizes for potential skipping events, derived from input annotation file, are indicated by the dashed ‘Annotated’ line.

Figure 6.

SpliceSiteScoring. Average splice site scores for upstream exon donor, skipped exon acceptor, skipped exon donor and downstream exon acceptor were determined for events with increased exon skipping (FDR < 0.0005) in RBP knockdowns and indisulam or ms023 treated cells. Only experiments with greater than 100 significant events were considered. A) Ratio of splice site scores for significant increased exon skipping events (–IncDiff) to significant decreased exon skipping events (+IncDiff). B) Average splice site scores for increased skipping events. Dashed lines provide reference scores for potential skipping configurations derived from an input annotation file.

Figure 7.

SpliceCompare. Comparative analysis of common increased exon skipping events (FDR < 0.0005) for RBP knockdowns and pharmacologic splicing inhibition.

Figure 8.

Neopeptide prediction from SETranslateNMD. Numbers of neopeptides longer than 12 amino acids long predicted from frameshifts generated by increased exon skipping (FDR < 0.0005) across the top 23 RBP knockdowns and pharmacologic splicing inhibition.

Figure 9.

NMD prediction from SETranslateNMD. Fraction predicted NMD isoforms are plotted for NMD inhibition and RBP knockdown experiments with a minimum of 100 increased and 100 decreased (FDR < 0.0005) skipped exon isoform level changes. SE = skipped exon.

Figure 10.

Predicting the functional impact of exon skipping. Pathway analysis of genes with increased exon skipping predicted to undergo NMD (upper panel) and with skipped conserved domain sequences (middle panel) were analyzed by Enrichr (https://maayanlab.cloud/Enrichr/) with the BioPlanet 2019 pathway database. Gene set enrichment analysis (GSEA) of gene expression (lower panel) was performed using E2F and MYC target signatures.