A Comprehensive Guide to Alternative Splicing Analysis with Bulk RNA-Seq: From Foundational Concepts to Clinical Applications

Abstract

This article provides a complete roadmap for researchers and drug development professionals conducting alternative splicing (AS) analysis using bulk RNA sequencing. It covers the foundational biology of AS and its significance in disease and development, explores the landscape of computational tools and methodologies for detecting splicing variations, addresses critical troubleshooting and optimization strategies for robust experimental design, and outlines best practices for validation and comparative analysis. By integrating the latest methodological advancements and empirical findings, this guide aims to empower scientists to generate biologically accurate and clinically relevant insights from splicing data, ultimately advancing drug discovery and precision medicine.

Understanding Alternative Splicing: Biological Significance and Transcriptome Diversity

The central dogma of molecular biology outlines the flow of genetic information from DNA to RNA to protein. A critical, regulated juncture in this pathway is the processing of precursor messenger RNA (pre-mRNA) into mature mRNA, a process dominated by splicing. In eukaryotes, the majority of multi-exon genes do not produce a single mRNA transcript but undergo alternative splicing (AS), where different combinations of exons are joined to generate multiple distinct isoforms from a single gene locus [1]. This process, coupled with mechanisms like alternative polyadenylation (APA), dramatically expands the informational capacity of the genome, serving as a major driver of proteome diversity and functional complexity [2].

Understanding the principles governing AS is not merely an academic pursuit but a foundational element of modern molecular biology and therapeutic development. In humans, over 95% of multi-exon genes are subject to AS, and dysregulation of splicing is a hallmark of numerous diseases, including cancers and neurodegenerative disorders [3] [4]. For researchers and drug development professionals, analyzing AS patterns provides critical insights into cellular states, disease mechanisms, and potential therapeutic targets. This document frames the core principles and methodologies for AS analysis within the context of bulk RNA sequencing (RNA-seq) research, providing a bridge from the fundamental mechanics of spliceosome action to the functional interpretation of proteomic outcomes.

Foundational Principles of Pre-mRNA Splicing and Regulation

2.1 The Spliceosome and Splicing Signals Pre-mRNA splicing is catalyzed by the spliceosome, a dynamic mega-dalton complex composed of five small nuclear ribonucleoproteins (snRNPs) and numerous associated proteins [1]. Its function is to recognize specific conserved sequence motifs within the pre-mRNA:

• The 5′ splice site (donor site)
• The 3′ splice site (acceptor site)
• The branch point sequence
• The polypyrimidine tract [1] Through a series of coordinated assembly steps and transesterification reactions, the spliceosome excises introns and ligates adjacent exons to form the mature mRNA.

2.2 Modes of Alternative Splicing AS subverts constitutive splicing patterns to produce variant transcripts. The major canonical types include:

• Cassette Exon (Exon Skipping): The most common and highly conserved type in mammals, where an exon is either included or skipped [3].
• Alternative 5′ or 3′ Splice Sites: Selection of different donor or acceptor sites, lengthening or shortening an exon.
• Intron Retention (IR): The failure to remove an intron, which is particularly prevalent in plants and linked to stress responses, and also significant in mammalian systems like cell fate transitions [1] [5].
• Mutually Exclusive Exons: Selection of one exon from a pair or set, often encoding structurally distinct protein domains.
• Alternative First or Last Exons: Changes in transcription start or polyadenylation sites, altering the N- or C-terminus of the protein [2].

2.3 Cis and Trans Regulation of Splicing Choices The selection of splice sites is not automatic but is tightly regulated by a combination of cis-acting RNA sequence elements and trans-acting factors.

• Cis-Regulatory Elements: These include exonic splicing enhancers (ESEs), exonic splicing silencers (ESSs), and their intronic counterparts (ISEs and ISSs). They are short, degenerate sequences that serve as binding platforms for regulatory proteins [1].
• Trans-Acting Factors: The primary regulators are RNA-binding proteins (RBPs), most notably the Serine/Arginine-rich (SR) proteins and heterogeneous nuclear ribonucleoproteins (hnRNPs). SR proteins typically bind enhancers and promote exon inclusion by recruiting core spliceosomal components, while hnRNPs often bind silencers and promote exon exclusion [1] [3]. The dynamic balance and tissue-specific expression of these factors ultimately determine the splicing outcome, integrating signals from development, cellular stress, and disease states.

Diagram: Core Splicing Regulation Pathway

Technological and Computational Frameworks for Analysis

3.1 Sequencing Platforms for Splicing Detection Bulk RNA-seq is the cornerstone for transcriptome-wide AS analysis. The choice between short-read and long-read technologies presents a trade-off between scale, accuracy, and isoform resolution.

Table 1: Comparison of RNA Sequencing Platforms for AS Analysis [1] [6]

Feature	Short-Read (e.g., Illumina)	Long-Read (PacBio HiFi)	Long-Read (Oxford Nanopore)
Read Length	Short (50-300 bp)	Long (1-10+ kb)	Very Long (1-100+ kb)
Primary Template	cDNA	cDNA	cDNA or direct RNA
Base Accuracy	Very High (>99.9%)	Very High (HiFi >99.9%)	Moderate (~96-99%)
Key Strength	High quantitative power, low cost per base	Full-length isoform resolution with high accuracy	Direct RNA seq, detects modifications, extreme length
Key Limitation	Indirect isoform inference; misses complex loci	Higher cost, lower throughput than short-read	Higher error rate can complicate splice site mapping
Best for AS	Quantifying known AS events (Psi) across many samples	De novo isoform discovery, complex loci	Detecting RNA modifications affecting splicing, novel isoforms

3.2 Computational Tools for AS Identification and Quantification Analyzing RNA-seq data to detect AS requires specialized computational tools that calculate metrics like Percent Spliced-In (PSI or Ψ), which quantifies the inclusion level of a particular exon or splice junction.

Table 2: Selected Computational Methods for AS Analysis from Bulk RNA-seq Data

Tool Name	Primary Function	AS Events Detected	Key Principle / Approach	Ref.
rMATS	Differential AS analysis	SE, MXE, A5SS, A3SS, RI	Uses a hierarchical model on junction-spanning and exon reads to compare PSI between groups.	[2]
SUPPA2	Differential AS & APA	Multiple AS events & APA	Calculates PSI from transcript expression; fast pattern generation across conditions.	[2] [7]
LeafCutter	Differential intron splicing	Intron clusters (splicing changes)	Identifies variable intron excision events without reliance on annotated splice sites.	[2]
MAJIQ	Quantification & visualization	Complex local splicing variations	Builds splicing graphs to model and quantify confidence in PSI changes.	[2]
Whippet	Differential AS analysis	Multiple AS events	Uses a lightweight probabilistic model leveraging gene annotation.	[2]

3.3 Experimental Validation Protocols Findings from computational analysis must be validated with targeted molecular biology assays.

• Protocol 1: Semi-Quantitative RT-PCR for Isoform Detection
  • Primer Design: Design primers in constitutive exons flanking the alternative region to amplify all isoforms.
  • cDNA Synthesis: Use high-quality, DNase-treated RNA and reverse transcriptase with random hexamers or oligo-dT.
  • PCR Amplification: Use a high-fidelity polymerase. Perform cycle titration to ensure amplification remains in the exponential phase for semi-quantitative comparison [1].
  • Analysis: Separate products via agarose or polyacrylamide gel electrophoresis. Bands can be excised and sequenced for confirmation. For higher resolution, use capillary electrophoresis (e.g., Agilent Bioanalyzer) [1].
• Protocol 2: Digital Droplet PCR (ddPCR) for Absolute Isoform Quantification
  • Assay Design: Design TaqMan probes or EvaGreen assays specific to each isoform (e.g., spanning an exon-exon junction unique to that isoform).
  • Partitioning & PCR: The reaction mix is partitioned into ~20,000 nanoliter-sized oil droplets. Each droplet acts as an independent PCR reactor [1].
  • End-Point Reading & Quantification: Droplets are read as positive or negative for fluorescence. Absolute copy numbers are calculated using Poisson statistics, without the need for a standard curve, offering high precision for low-abundance isoforms [1].
• Protocol 3: Minigene Splicing Reporter Assay
  • Construct Cloning: Clone the genomic region of interest (containing the alternative exon and its flanking introns) into an exon-trapping vector between two constitutive exons.
  • Transfection: Transfect the construct into relevant cell lines.
  • RNA Isolation & Analysis: Isolve RNA 24-48 hours post-transfection and analyze splicing patterns via RT-PCR specific to the minigene transcript. This assay isolates the splicing regulation of a specific locus from genomic context [1].

Diagram: Bulk RNA-seq AS Analysis Workflow

From Splicing to Functional Proteome Consequences

4.1 Predicting and Modeling Structural Impacts The ultimate consequence of AS is a change in protein sequence, which can have profound effects on structure and function. Computational protein structure prediction tools like AlphaFold2 now enable systematic investigation of these impacts. Studies predicting structures for thousands of human splice isoforms reveal that AS can alter secondary structure composition, surface charge distribution, and protein compactness (radius of gyration) [4]. Crucially, AS often buries or exposes post-translational modification (PTM) sites, fundamentally altering regulatory potential. For example, alternative splicing of the BAX gene modulates exposure of phosphorylation sites, affecting its pro-apoptotic activity [4].

4.2 Functional Enrichment and Biological Interpretation To move from a list of differentially spliced genes to biological insight, functional enrichment analysis is essential. Genes undergoing AS during specific processes—like the fast chemical reprogramming (FCR) of somatic cells—are enriched in stage-specific pathways. In late FCR stages, spliced genes relate to ubiquitin-mediated proteolysis, while in established induced pluripotent stem cells (iPSCs), they are enriched in signaling pathways like mTOR and PI3K-Akt [5]. This indicates AS plays stage-specific regulatory roles. Similarly, cross-species analysis has identified AS events in genes related to stress response and neuronal function as being associated with maximum lifespan (MLS), suggesting a role for splicing in longevity [3].

4.3 Integration with Single-Cell and Spatial Context While bulk RNA-seq measures average splicing patterns, emerging technologies contextualize isoform expression. Long-read single-cell RNA-seq (e.g., Nanopore) allows isoform resolution at the cellular level, revealing cell-type-specific splicing and intra-cell isoform heterogeneity [8]. Integrating bulk-derived AS insights with single-cell atlases (e.g., Tabula Sapiens) can predict in which cell types a structurally distinct isoform is expressed, linking molecular form to cellular function [4].

Table 3: Summary of Key Experimental Validation Methods

Method	Primary Application	Key Metric	Advantages	Limitations
RT-PCR / Capillary Electrophoresis	Detect & semi-quantify known isoforms.	Amplicon size / peak area.	Low-cost, rapid, high-specificity.	Low throughput; semi-quantitative.
Digital Droplet PCR (ddPCR)	Absolute quantification of specific isoforms.	Copies/μL of target isoform.	High precision, absolute quantification, no standard curve needed.	Low multiplexing; requires specific assay design.
Minigene / Splicing Reporter	Functional testing of cis-elements and trans-factors.	PSI from reporter transcript.	Isolates regulatory sequence; can test mutants.	May lack full genomic/chromatin context.

Diagram: From Splicing Event to Functional Hypothesis

Table 4: Key Research Reagent Solutions for Splicing Analysis

Item	Function / Application	Key Considerations
High-Quality RNA Isolation Kits	Obtain intact, degradation-free RNA for accurate isoform representation.	Prioritize methods with effective DNase treatment and inhibitors of RNases.
Ribonuclease H (RNase H)-deficient Reverse Transcriptases	Generate full-length cDNA without degrading the RNA template during first-strand synthesis.	Essential for long amplicons in validation assays.
High-Fidelity DNA Polymerases	Amplify cDNA for validation (RT-PCR) with minimal errors, ensuring accurate isoform sequence.	Critical for downstream sequencing of amplicons.
Spliceosome & RBP Inhibitors	Mechanistic studies (e.g., Pladienolide B for SF3B1).	Tools to perturb splicing and study immediate downstream effects.
Isoform-Specific Antibodies	Validate protein-level consequences of AS via western blot or immunofluorescence.	Rare; often require custom generation against isoform-unique peptides.
Long-Read Sequencing Kits	For full-length isoform discovery (PacBio) or direct RNA modification analysis (Nanopore).	Choice depends on need for accuracy vs. read length/modification detection [6].
Bioinformatics Pipelines & Software	Align reads, quantify expression, identify and differentially test AS events (e.g., rMATS, SUPPA2).	Containerized versions (Docker/Singularity) ensure reproducibility [2].
Splicing-Focused Databases	Query known isoforms, conservation, and tissue-specific patterns (e.g., VastDB, Ensembl).	Provide context for novel findings [5].

Alternative splicing (AS) is a fundamental post-transcriptional process that enables a single gene to produce multiple mRNA isoforms by differentially including or excluding exonic and intronic regions [9]. This mechanism is a major source of proteomic diversity and a critical regulator of development, tissue specificity, and cellular differentiation [10]. In humans, over 95% of multi-exon genes undergo alternative splicing [11], and its dysregulation is implicated in numerous diseases, including cancer, neurodegenerative disorders, and muscular dystrophies [12] [13].

The advent of bulk RNA sequencing (RNA-seq) has revolutionized the study of AS, providing a high-throughput, quantitative platform for discovering and quantifying splicing variations across biological conditions. Within the framework of a bulk RNA-seq research thesis, precise identification and quantification of splicing events—primarily exon skipping (ES), intron retention (IR), and alternative 5′/3′ splice site (ASSS) selection—are paramount. These events represent distinct mechanistic outcomes with unique biological implications and technical challenges for detection. This application note details the molecular basis, detection methodologies, and analytical protocols for these major splicing types, providing a practical guide for researchers and therapeutic developers.

Molecular Mechanisms and Biological Significance

Splicing is catalyzed by the spliceosome, a dynamic complex of small nuclear ribonucleoproteins (snRNPs) and associated proteins that recognize conserved sequences at exon-intron boundaries: the 5′ splice site (SS), branch point sequence (BPS), polypyrimidine tract (PPT), and 3′ SS [10] [9]. The selection of these sites is regulated by cis-acting RNA elements (enhancers and silencers) and trans-acting RNA-binding proteins (e.g., SR proteins and hnRNPs), which together form a complex "splicing code" [11] [9].

Exon Skipping (ES), or cassette exon, is the most prevalent AS event in mammals [11] [9]. It involves the complete exclusion of an exon from the mature mRNA, along with its flanking introns. ES can lead to the loss of protein domains, affect protein-protein interactions, or cause a frameshift. It is a primary therapeutic target; for example, antisense oligonucleotides (AONs) are designed to skip exons harboring deleterious mutations in DMD to restore the reading frame in Duchenne Muscular Dystrophy (DMD) [12].

Intron Retention (IR) occurs when an intron is not removed and remains in the mature mRNA. It is the most common AS event in plants but less frequent in mammals, where it is often associated with nonsense-mediated decay (NMD) or regulation of gene expression [11] [9]. Retained introns can introduce premature termination codons or alter the protein's C-terminus if located in the final exon.

Alternative 5′ or 3′ Splice Site (ASSS) Selection involves the use of different donor (5′) or acceptor (3′) sites within an exon or intron. This leads to the extension or truncation of an exon, subtly altering the protein's coding sequence [11] [13]. The selection between competing splice sites is influenced by their intrinsic strength, the local concentration of splicing factors, and a competitive mechanism where nearby sites influence each other's usage [14].

The following table summarizes the key characteristics of these major event types.

Table 1: Characteristics of Major Alternative Splicing Events

Event Type	Description	Prevalence in Humans	Potential Functional Impact	Therapeutic Example
Exon Skipping (ES)	Complete exclusion of an exon from the mature transcript.	~30% of AS events; most common pattern [11].	Domain loss, altered protein function, frameshift.	AONs (Eteplirsen) to skip exon 51 in DMD for DMD [12].
Intron Retention (IR)	An intron remains in the final mRNA.	Less common; frequent in untranslated regions (UTRs) [11].	Introduction of PTCs (leading to NMD), altered C-terminus, regulatory non-coding RNAs.	Targeted for degradation or modulation of gene expression levels.
Alt. 5′/3′ Splice Site (ASSS)	Use of alternative donor (5′) or acceptor (3′) sites.	~25% of AS events [11].	Subtle insertion/deletion in coding sequence, minor protein domain alteration.	Correction of cryptic splice site usage caused by deep-intronic mutations [13].

Detection and Quantification from Bulk RNA-Seq Data

Accurate analysis of AS from short-read bulk RNA-seq data involves specific computational strategies. The core metric for quantification is the Percent Spliced In (PSI or Ψ), which estimates the relative inclusion level of an exon, intron, or alternative splice site [15].

3.1. Computational Approaches and Tools Analysis pipelines can be categorized into splice-junction-centric and exon-centric methods. Junction-centric tools (e.g., rMATS, MAJIQ, LeafCutter) quantify reads spanning splice junctions to infer PSI. Exon-centric tools (e.g., DEXSeq) model reads covering exonic regions to test for differential usage. A significant advancement is the incorporation of exon-exon junction reads alongside exon counts, which resolves double-counting issues and enhances power, particularly for detecting ASS and IR events [16].

For large-scale or heterogeneous datasets (e.g., GTEx, TCGA), tools like MAJIQ v2 are essential. It employs a Local Splicing Variation (LSV) framework that captures complex events and unannotated junctions, and uses a Bayesian model to quantify PSI and confidence intervals. Its "HET" module applies non-parametric tests to handle sample heterogeneity effectively [15].

3.2. Performance of Detection Methods The choice of tool significantly impacts detection capability. Benchmarking studies show that methods incorporating junction information (DEJU) outperform traditional exon-only (DEU) approaches [16].

Table 2: Performance Comparison of Splicing Detection Methods for Different Event Types [16]

Splicing Event Type	Key Challenge in Detection	Recommended Method (High Power & Controlled FDR)	Method with Notable Limitations
Exon Skipping (ES)	Distinguishing skipping from constitutive exon.	DEJU-edgeR/limma, rMATS, MAJIQ v2.	Basic DEXSeq (lower power for complex genes).
Intron Retention (IR)	Differentiating true retention from pre-mRNA contamination or reads from unprocessed RNA.	MAJIQ v2, DEJU-edgeR/limma (specifically designed for IR).	Most standard DEU workflows (fail to detect IR) [17].
Alt. Splice Site (ASSS)	Identifying subtle changes in exon boundaries.	DEJU-edgeR/limma, MAJIQ v2 (LSV framework).	Exon-only DEU approaches (lack junction resolution).

3.3. Visualization of Splicing Analysis Workflow The following diagram illustrates a generalized computational workflow for differential splicing analysis from raw RNA-seq reads to biological insight, integrating steps from alignment, quantification, statistical testing, and visualization.

Diagram Title: Computational Workflow for Differential Splicing Analysis

Detailed Experimental and Analytical Protocols

4.1. Protocol: Validating Exon Skipping Using RT-PCR and Sequencing This protocol is used to experimentally confirm in silico predictions of exon skipping, common in therapeutic development for diseases like DMD [12].

• RNA Isolation & QC: Extract total RNA from treated and control cells/tissues using TRIzol. Assess integrity (RIN > 8.0) via Bioanalyzer.
• cDNA Synthesis: Perform reverse transcription using 1 µg of total RNA, oligo(dT) primers, and a reverse transcriptase.
• PCR Amplification: Design primers in the exons flanking the predicted skipped exon. Perform PCR with a high-fidelity polymerase. Include a control amplifying a constitutively spliced gene.
• Gel Electrophoresis: Resolve PCR products on a 2-3% agarose gel. The skipped isoform will produce a smaller band than the full-length product.
• Purification & Sequencing: Excise gel bands, purify the DNA, and perform Sanger sequencing to confirm the exact junction.

4.2. Protocol: Differential Splicing Analysis with MAJIQ v2 on a Large Dataset This protocol is designed for robust analysis of heterogeneous data across many samples, such as different brain regions or cancer subtypes [15].

• Data Preparation: Compile BAM alignment files for all samples. Prepare a sample table defining group affiliations (e.g., condition, tissue).
• Build Splice Graphs: Run the MAJIQ builder to create a unified splice graph for each gene, incorporating annotated and de novo junctions.

• Quantify Splicing Variations: Run the MAJIQ quantifier to compute PSI (Ψ) distributions for all LSVs in each sample group.

• Identify Differential LSVs: Apply the MAJIQ HET test to identify LSVs with significant ΔPSI between groups (e.g., |ΔΨ| > 0.2, probability > 0.95).
• Visualize Results: Use VOILA v2 to visualize complex splicing events, generating splice graphs and posterior distribution plots for top hits.

4.3. Protocol: Targeted Analysis of Alternative Splice Sites via Amplicon-Seq This protocol is useful for deep investigation of specific ASS events predicted by tools like those in [14] or for screening patient variants.

• Design Amplification Strategy: Design long-range PCR primers to amplify the genomic region encompassing the alternative exon and its flanking constitutive exons.
• Library Preparation & Sequencing: Amplify the region from cDNA, prepare an Illumina amplicon library, and perform paired-end 2x150 bp or 2x250 bp sequencing for sufficient junction overlap.
• Mapping & Junction Analysis: Map reads to the reference genome using a sensitive aligner (e.g., BWA-MEM). Extract and count all unique splice junction combinations.
• Quantify Site Usage: Calculate the PSI for each competing 5' or 3' splice site by dividing its junction read count by the sum of all junction reads for that alternative splicing event.

The Scientist's Toolkit: Research Reagent Solutions

Table 3: Essential Reagents and Tools for Splicing Research

Item	Function / Description	Example Use Case / Note
Splice-Aware Aligner	Aligns RNA-seq reads across splice junctions.	STAR [16] or HISAT2 for generating BAM files for downstream analysis.
Splicing Quantification Software	Quantifies PSI and detects differential splicing.	MAJIQ v2 [15] for complex/heterogeneous data; rMATS for standard case-control studies.
Antisense Oligonucleotides (AONs)	Chemically modified RNAs/DNAs that bind pre-mRNA to modulate splicing.	Induce exon skipping (Eteplirsen for DMD) [12] or block cryptic splice sites.
U7 snRNA Vector System	Engineered snRNA for sustained in vivo expression of antisense sequences.	Used in recombinant AAV vectors for long-term exon skipping therapy in muscle [12].
Splicing Reporter Minigenes	Plasmid constructs containing a genomic region of interest cloned into an exon-intron-exon reporter.	Validate the impact of genetic variants or ASOs on splicing patterns in vitro [13].
RT-PCR Reagents & Gel Electrophoresis	Standard molecular biology tools for isoform detection.	Rapid, low-cost validation of major splicing changes (e.g., exon skipping) [12].

The precise analysis of exon skipping, intron retention, and alternative splice site usage is a cornerstone of modern transcriptomics and a critical component of a bulk RNA-seq research thesis. The integration of sophisticated computational tools like MAJIQ v2 [15] and DEJU [16] with targeted experimental validation provides a robust framework for discovering disease-relevant splicing dysregulation.

The translational impact is profound. Exon skipping therapies for DMD exemplify successful "splice-switching," with approved AONs and advanced AAV-U7snRNA vectors in clinical trials [12]. Furthermore, the systematic annotation of splice-disruptive variants is enhancing diagnostic yield and revealing new targets for RNA-targeted therapies in neurological and other disorders [13]. As long-read sequencing and single-cell technologies mature, they will further refine our understanding of splicing diversity. However, bulk RNA-seq remains essential for its cost-effectiveness, statistical power in cohort studies, and well-established analytical pipelines, continuing to be an indispensable tool for researchers and drug developers aiming to decode and manipulate the splicing code for therapeutic benefit.

The Critical Role of Splicing in Development, Cell Differentiation, and Disease

Biological Context and Analytical Imperative

Alternative splicing (AS) is a fundamental post-transcriptional regulatory mechanism, affecting over 90% of multi-exon human genes and enabling a limited genome to produce a vast diversity of protein isoforms [8]. In development, precise splicing patterns dictate cell fate determination and lineage specification. During differentiation, splicing shifts remodel the cellular proteome to support specialized functions. Dysregulation of this precise control is a hallmark of disease; aberrant splicing contributes directly to pathologies including cancer, intellectual disability, and autism spectrum disorders [8].

Bulk RNA sequencing (RNA-seq) remains a cornerstone for profiling these splicing dynamics across tissues and patient cohorts. It provides a population-average view essential for identifying consistent splicing signatures associated with biological states or clinical outcomes. The core analytical challenge lies in accurately quantifying transcript isoforms from short sequencing reads that may ambiguously map to multiple splice variants [18]. Moving from gene-level expression to isoform-resolution analysis is therefore critical for a mechanistic understanding of how splicing influences biology and disease.

Foundational Analytical Framework for Splicing from Bulk RNA-seq

A robust analysis of alternative splicing from bulk RNA-seq data requires a multi-step computational workflow designed to handle uncertainty in read assignment and quantify splice variant abundance.

Core Workflow Steps:

• Read Alignment & Quantification: Sequencing reads are aligned to a reference genome using a splice-aware aligner (e.g., STAR). Subsequently, tools like Salmon or RSEM perform transcript-level quantification, using statistical models to resolve the uncertainty when reads map to multiple isoforms [18].
• Splicing Event Identification: Quantified transcripts are analyzed to detect specific types of splicing events (e.g., exon skipping, alternative 5'/3' splice sites).
• Differential Splicing Analysis: Statistical testing is applied to identify events whose usage (e.g., Percent Spliced In, PSI) changes significantly between experimental conditions or patient groups.

The following workflow diagram outlines this generalized pipeline for differential splicing analysis.

Diagram: Bulk RNA-seq Splicing Analysis Workflow.

Experimental Protocol: Splicing-Focused Bulk RNA-seq Library Preparation and Sequencing

This protocol details the steps for generating sequencing libraries suitable for alternative splicing analysis from total RNA [19] [18].

Materials and Equipment

• Total RNA samples (RIN > 8.0 recommended).
• Poly(A) mRNA Magnetic Isolation Kit (e.g., NEBNext) [19].
• RNA-seq Library Prep Kit for Illumina (e.g., NEBNext Ultra II).
• SPRIselect Beads or equivalent for size selection and clean-up.
• Thermal cycler, magnetic stand, and microcentrifuge.
• High-sensitivity DNA assay (e.g., Agilent Bioanalyzer).
• Illumina sequencing platform (e.g., NextSeq, NovaSeq).

Step-by-Step Procedure

A. mRNA Enrichment and Fragmentation

• Isolate polyadenylated mRNA from 50 ng – 1 µg of total RNA using poly(A) magnetic beads. This step depletes ribosomal RNA.
• Elute the purified mRNA in nuclease-free water.
• Fragment the mRNA using divalent cations under elevated temperature (e.g., 94°C for 5-15 minutes) to generate fragments of 200-300 base pairs.

B. First and Second Strand cDNA Synthesis

• Synthesize first-strand cDNA using random hexamer priming and reverse transcriptase.
• Synthesize the second-strand cDNA using DNA Polymerase I and RNase H, creating double-stranded cDNA.

C. Library Construction

• End Repair/A-Tailing: Repair fragment ends and add a single 'A' nucleotide to the 3' ends to facilitate adapter ligation.
• Adapter Ligation: Ligate indexed sequencing adapters with a single 'T' overhang to the cDNA fragments.
• Library Amplification: Perform 8-12 cycles of PCR to enrich for adapter-ligated fragments and incorporate full-length adapter sequences.

D. Quality Control and Pooling

• Clean the amplified library using SPRIselect beads at a 0.8x-1.0x bead-to-sample ratio to remove short fragments and reactants.
• Quantify the final library using a fluorometric method (e.g., Qubit) and assess the size distribution via Bioanalyzer. Expect a peak at ~300-400 bp.
• Pool equimolar amounts of uniquely indexed libraries for multiplexed sequencing.

E. Sequencing

• Sequence the pooled library on an Illumina platform. For optimal splicing analysis, use paired-end sequencing (e.g., 2x150 bp) to generate reads long enough to span splice junctions [18].

Computational Tools for Splicing Analysis: A Comparative Guide

A wide array of computational tools has been developed to detect and quantify alternative splicing events from RNA-seq data. They can be broadly categorized by their primary function and methodological approach [2].

Table 1: Computational Tools for Alternative Splicing Analysis from Bulk RNA-seq Data

Tool Name	Primary Function	Method/Approach	Key Strength	Reference
rMATS	Differential AS detection	Uses a hierarchical model to compare splice junction counts.	Models biological variance; handles replicates well.	[2]
SUPPA2	Differential AS & APA	Calculates PSI from transcript expression; fast profiling.	High speed; integrates alternative polyadenylation.	[2]
LeafCutter	Differential intron splicing	Clusters intron excision events without prior annotation.	Detects novel and complex splicing events.	[2]
MAJIQ	Quantification & differential AS	Builds local splice graphs to model complex splicing variations.	Resolves complex splicing variations effectively.	[2]
DEXSeq	Differential exon usage	Tests for differential usage of exonic bins between conditions.	Generalizable to any genomic feature.	[2]

The analysis can target seven major types of alternative splicing events, each with distinct functional consequences for the resulting mRNA and protein.

Diagram: Major Types of Alternative Splicing Events.

Successful splicing analysis relies on both wet-lab reagents and bioinformatic resources. The following table details key solutions.

Table 2: Research Reagent Solutions for Splicing Analysis

Category	Item / Solution	Function & Importance	Example / Note
Wet-Lab	Poly(A) mRNA Selection Beads	Enriches for mature, polyadenylated mRNA, removing rRNA that would dominate sequencing.	NEBNext Poly(A) mRNA Magnetic Beads [19].
Wet-Lab	Stranded Library Prep Kit	Preserves strand-of-origin information, crucial for accurate transcript annotation and splicing detection.	Illumina Stranded mRNA Prep.
Wet-Lab	High-Fidelity Polymerase	Reduces PCR errors during library amplification, ensuring sequence fidelity for variant detection.	Q5 Hot Start Polymerase.
Bioinformatic	Splice-Aware Aligner	Aligns RNA-seq reads across splice junctions by allowing gaps in the alignment.	STAR, HISAT2 [18].
Bioinformatic	Transcript Quantifier	Estimates isoform abundance from aligned reads, resolving multimapping ambiguity.	Salmon (alignment-based mode), RSEM [18].
Bioinformatic	Splicing Analysis Tool	Performs statistical testing for differential splicing between sample groups.	rMATS, SUPPA2 (see Table 1) [2].
Reference	High-Quality Annotation	Provides a reference set of known transcripts and splice junctions for quantification.	GENCODE, Ensembl GTF files.

From Analysis to Therapeutic Insight

Integrating splicing analysis into a broader research thesis provides a direct path to mechanistic and translational discovery. In cancer research, identifying tumor-specific splice variants (neoantigens) can reveal novel therapeutic targets and immunotherapies. For neurodevelopmental disorders, elucidating splicing dysregulation in genes like RBFOX or MECP2 can pinpoint critical windows for intervention.

The emerging integration of single-cell and spatial long-read sequencing, as exemplified by tools like Longcell which addresses UMI recovery and isoform quantification challenges [8], promises to map splicing heterogeneity within tissues. This refines our understanding from population averages to cellular microenvironments, revealing how splicing contributes to cell fate decisions in development or tumor microevolution in disease. Ultimately, the systematic application of these bulk RNA-seq protocols and analytical frameworks is foundational for dissecting the critical role of splicing across biology and advancing RNA-targeted therapeutics.

The Central Role of Alternative Splicing in Human Disease and Therapy

Alternative splicing (AS) is a fundamental post-transcriptional mechanism through which a single gene can generate multiple mRNA isoforms, leading to proteomic diversity. Recent analyses indicate that approximately 95% of human multi-exon genes undergo alternative splicing [20]. This process is not merely a source of diversity; its dysregulation is a direct contributor to disease pathogenesis. It is estimated that 15-50% of human genetic diseases involve mutations that disrupt normal splicing patterns, affecting everything from core splice sites to splicing regulatory elements (SREs) [20] [21]. Consequently, the precise patterns of splicing in a tissue or cell type are not just molecular noise but are rich sources of disease-specific targets and clinically actionable biomarkers.

In cancer, for instance, splicing dysregulation can create isoforms that drive tumorigenesis, promote metastasis, or confer drug resistance. The splicing factor/binding site axis represents a promising but complex therapeutic frontier. As evidenced by the development of splicing modulator compounds (SMCs) like BPN-15477, targeted correction of pathogenic splicing is a viable therapeutic strategy [22]. This compound was shown to correct aberrant splicing in the ELP1 gene linked to familial dysautonomia and was predicted via a deep learning model to be capable of correcting splicing defects in 155 other disease genes [22]. For drug discovery professionals, this underscores a dual opportunity: splicing variants can be the therapeutic target themselves, or their unique expression profiles can serve as biomarkers for patient stratification, treatment response, and disease monitoring.

Computational Landscape for Splicing Analysis from Bulk RNA-Seq

The advent of high-throughput RNA sequencing (RNA-seq) has been transformative for splicing analysis. Bulk RNA-seq analysis of splicing has evolved from quantifying known isoforms to detecting complex, unannotated local splicing variations (LSVs) across heterogeneous sample sets. The core metric for quantification is the Percent Spliced In (PSI or Ψ) value, which measures the relative inclusion of an exon or junction [15].

A critical challenge in the field is the accurate analysis of large-scale, heterogeneous datasets (e.g., from consortia like GTEx). Traditional tools that assume homogeneous sample groups often fail or lose power with such data. Next-generation pipelines like MAJIQ v2 address this by introducing heterogeneous (HET) test statistics that quantify PSI per sample before group comparison, enhancing robustness and power [15]. Furthermore, tools must grapple with transcriptome complexity, where a significant fraction of splicing variations involve unannotated junctions or complex patterns not captured by simple event types [15] [23].

Table 1: Common Types of Alternative Splicing Events and Their Features

Splicing Event Type	Description	Approximate Frequency in Humans [20]	Potential Functional Impact
Exon Skipping (ES)	An exon is either included or skipped in the mature mRNA.	~38.4% (Most common)	Can lead to loss or gain of protein domains, altering function or stability.
Alternative 3' Splice Site (A3SS)	Selection between different splice sites at the 3' end of an upstream exon.	~18.4%	Changes the C-terminal boundary of an upstream exon, potentially affecting protein sequence.
Alternative 5' Splice Site (A5SS)	Selection between different splice sites at the 5' end of a downstream exon.	~7.9%	Changes the N-terminal boundary of a downstream exon.
Intron Retention (IR)	An intron is retained in the mature transcript.	~2.8%	Often introduces a premature termination codon (PTC), leading to nonsense-mediated decay (NMD) or a protein with an inserted sequence.
Mutually Exclusive Exons (MXE)	One of two adjacent exons is included, but never both.	Part of ~32.4% "Other"	Provides a switch between two distinct protein modules.

Table 2: Selected Computational Tools for Bulk RNA-seq Splicing Analysis

Tool Name	Core Methodology	Key Strength	Best Suited For
MAJIQ v2 [15]	Bayesian modeling of Local Splicing Variations (LSVs); HET non-parametric tests.	Handles large, heterogeneous datasets; detects complex & de novo events.	Large cohort studies (e.g., GTEx), differential splicing across many conditions.
rMATS [24]	Statistical modeling of reads spanning junction pairs for defined event types.	High precision for classical event types; widely used and validated.	Standard differential splicing analysis between controlled experimental groups.
SpliceSeq [25]	Alignment to splice graphs; functional consequence prediction.	Intuitive visualization; links splicing changes to protein function via UniProt.	Prioritizing splicing events with likely functional, pathogenic impact.
AStalavista [23]	Systematic identification and categorization of AS events from annotation.	Provides a universal "AS code" nomenclature; characterizes full splicing landscape.	Categorizing and comparing global splicing patterns across annotations or species.
IRFinder [24]	Specialized quantification of intron retention from RNA-seq.	High accuracy for detecting and quantifying IR events.	Studies where intron retention is a key focus (e.g., some cancers, viral infection).

Application Note 1: Protocol for Differential Splicing Analysis in Drug Response Studies

Objective: To identify splicing alterations associated with drug treatment or resistance using bulk RNA-seq data from cell lines or patient-derived samples.

Experimental Design:

• Sample Groups: Include treated vs. untreated controls, isogenic sensitive vs. resistant cell lines, or patient pre- vs. post-treatment samples. Minimum 3-4 biological replicates per group is strongly recommended.
• RNA-seq Library Preparation: Use poly-A selection to enrich for mature mRNA. Paired-end sequencing (e.g., 2x150 bp) with a depth of 50-100 million reads per sample is standard for splicing analysis to ensure sufficient junction coverage.
• Controls: Spike-in RNA controls (e.g., ERCC) can help monitor technical variance.

Bioinformatics Protocol (Using MAJIQ v2 as an example):

• Read Alignment and Preprocessing:
  • Align reads to the human reference genome (e.g., GRCh38) using a splice-aware aligner like STAR.
  • Process aligned BAM files: sort, index, and optionally remove duplicates.
• Build Splicing Graphs:
  • Run majiq build with a reference transcript annotation (e.g., GENCODE). This creates a splice graph for each gene, incorporating known and de novo detected junctions. For incremental analysis of large projects, MAJIQ v2 allows adding new samples without rebuilding from scratch [15].
• Quantify Splicing Variations:
  • Run majiq quant on your sample groups. This calculates the PSI (Ψ) and its posterior distribution for every Local Splicing Variation (LSV) in each sample/group. For heterogeneous data, specify the --het flag to use sample-wise quantification.
• Identify Differential Splicing (DS):
  • Run majiq deltapsi to compare groups (e.g., treated vs. control). This calculates ΔPSI (change in PSI) and the probability that |ΔPSI| > a threshold (e.g., 0.2 or 20%). Output is a list of LSVs with significant changes.
• Visualization and Interpretation:
  • Use VOILA v2 (the companion visualizer for MAJIQ) to inspect significant LSVs. It can generate simplified splice graphs showing PSI across groups for hundreds of samples [15].
  • Integrate gene ontology (GO) analysis on genes with significant DS events to identify enriched biological pathways affected by the drug.

Key Consideration: Always plan for experimental validation of top hits using an orthogonal method such as RT-PCR with capillary electrophoresis or nanostring nCounter on independent samples.

Application Note 2: Protocol for Identifying Splicing-Derived Biomarkers

Objective: To discover and prioritize tissue- or disease-specific splicing variants as potential diagnostic or prognostic biomarkers.

Workflow Overview:

• Discovery Cohort Analysis: Perform differential splicing analysis (as in Protocol 3.1) on a well-powered cohort comparing disease vs. normal tissues. Focus on events with large ΔPSI magnitude and high statistical significance.
• Prioritization Filter:
  • Specificity: Prioritize events that are tissue-specific or show minimal variation in normal tissues [20]. Tools like MAJIQ v2's "Modulizer" can help classify events [15].
  • Function: Use annotation (e.g., via SpliceSeq [25]) to prioritize events that alter critical protein domains, active sites, or structured motifs.
  • Detectability: Favor events where the splice junction boundaries are amenable to detection by qPCR, droplet digital PCR (ddPCR), or targeted RNA-seq assays.
• Biomarker Assay Development:
  • Design junction-spanning TaqMan assays or amplicon-seq probes specific to the biomarker isoform.
  • Optimize assay using a training set of samples.
• Validation: Test the biomarker assay on a large, independent validation cohort of patient samples, linking the splicing signature to clinical outcomes (e.g., survival, therapy response).

Diagram 1: Workflow for Splicing Biomarker Discovery & Validation

Table 3: Key Research Reagent Solutions for Splicing Studies

Category	Item / Resource	Function in Splicing Research	Example / Note
Wet-Lab Reagents	Poly-A Selection Beads	Isolate mature, poly-adenylated mRNA from total RNA for RNA-seq library prep.	Essential for standard RNA-seq protocols to focus on spliced transcripts.
	Ribonuclease H (RNase H)	Used in assays to specifically degrade RNA in DNA-RNA hybrids (R-loops), which can influence splicing.	Useful for studying R-loop mediated splicing regulation.
	Splicing Modulator Compounds (SMCs)	Small molecules that bind the spliceosome or splicing factors to influence splice site choice.	BPN-15477 [22]; tools for mechanistic studies & therapeutic leads.
Molecular Tools	Minigene Splicing Reporters	Plasmid-based systems containing a genomic region of interest to study splicing regulation in vivo.	Critical for validating the impact of mutations or SREs on splicing [21].
	Antisense Oligonucleotides (ASOs)	Chemically modified oligonucleotides that bind pre-mRNA to block or enhance splice site usage.	Used for functional validation and as a therapeutic modality (e.g., Nusinersen for SMA).
	Clontech SMARTER cDNA Kits	Generate high-quality cDNA from low-input RNA for isoform-specific PCR validation.	For validating RNA-seq findings via RT-PCR.
Computational Resources	SpliceAid-F / SpliceAid 2	Database of experimentally validated splicing factor binding motifs and sites [21].	Predicts SREs and interprets the effect of sequence variants.
	Ensembl / GENCODE Annotation	Comprehensive, regularly updated transcript annotation for the human genome.	Provides the reference "splice graph" for alignment and quantification.
	GTEx Portal	Repository of normal tissue RNA-seq data.	Serves as a critical baseline for assessing tissue-specificity of splicing events.

Advanced Frontiers: Integrating Splicing Analysis into Target Discovery

The future of splicing in drug discovery lies in integration and precision. Key frontiers include:

• Single-Cell Splicing Resolution: New computational frameworks like SCSES (2025) are overcoming the sparsity of scRNA-seq data to impute reliable PSI values, allowing the mapping of splicing heterogeneity within tumors and its correlation with drug-resistant cell states [24]. This can reveal rare, targetable subpopulations.
• AI-Driven Target Prediction: As demonstrated with BPN-15477, deep learning models (e.g., CNNs) trained on RNA-seq data from SMC-treated cells can identify responsive sequence signatures [22]. These models can then scan databases like ClinVar to predict which patient mutations are correctable by a given SMC, enabling a precision medicine approach.
• Integrative Multi-Omics: Combining splicing data with genetic variants (e.g., eQTLs for splicing - sQTLs), chromatin accessibility (ATAC-seq), and RBP binding (CLIP-seq) can elucidate the full regulatory network governing a pathogenic splicing event, revealing upstream therapeutic nodes.

Diagram 2: Therapeutic Targeting of a Splicing-Centric Disease Pathway

Concluding Synthesis

The systematic analysis of alternative splicing through bulk RNA-seq has matured from a basic research activity to a cornerstone of modern therapeutic discovery. The protocols and frameworks outlined here enable researchers to transition from identifying splicing correlates of disease to validating causal splicing targets and developing robust isoform-specific biomarkers. The integration of sophisticated computational tools (like MAJIQ v2 for bulk, SCSES for single-cell) with experimental validation creates a powerful pipeline for translating splicing biology into novel clinical assets. As the field progresses, the convergence of deep learning, single-cell genomics, and targeted splicing therapeutics promises to unlock a new class of precision medicines for a wide array of genetically defined disorders.

Bulk RNA-Seq as a Key Tool for Transcriptome-Wide Splicing Profiling

Bulk RNA sequencing (RNA-Seq) is a foundational technology for transcriptome-wide profiling of gene expression and, crucially, for the systematic analysis of alternative splicing (AS) [26]. By sequencing RNA from a population of cells, it provides an averaged but comprehensive snapshot of the transcriptome, enabling the identification and quantification of diverse RNA isoforms produced from a single gene locus [2]. Within a broader thesis on alternative splicing analysis, bulk RNA-Seq serves as the primary, high-throughput discovery tool. It allows researchers to detect splicing variations across different biological conditions—such as disease states, developmental stages, or drug treatments—on a genome-wide scale [15]. The ability to profile splicing in a transcriptome-wide manner is critical for understanding the complex regulatory networks underlying cellular differentiation, homeostasis, and disease pathogenesis, including cancer and autoimmune disorders [27] [28].

The analysis of alternative splicing via bulk RNA-Seq involves specific computational and statistical challenges distinct from differential gene expression analysis. The goal is to accurately quantify shifts in the relative abundance of mRNA isoforms. This is typically measured by metrics like the Percent Spliced In (PSI or Ψ), which represents the proportion of transcripts that include a particular exon or splice junction [15]. The power of bulk RNA-Seq in this context is its scalability and cost-effectiveness for profiling large sample cohorts, which is essential for achieving the statistical power needed to detect subtle but biologically significant splicing changes [29]. As part of an integrated analytical thesis, findings from bulk RNA-Seq often form the hypothesis-generating bedrock, which can be further validated with targeted assays or explored at single-cell resolution [27] [30].

Quantitative Foundations of Splicing Analysis

Catalog of Alternative Splicing Events

Alternative splicing generates transcriptomic diversity through several canonical patterns. The table below defines the primary AS event types analyzed in bulk RNA-Seq data and their functional implications [2].

Table 1: Primary Types of Alternative Splicing Events Detectable by Bulk RNA-Seq

Event Type	Description	Key Functional Consequence
Exon Skipping (Cassette Exon)	The complete inclusion or exclusion of an exon in the mature mRNA.	Most common in mammals; can lead to loss/gain of protein domains.
Alternative 5' Splice Site (Alt 5'SS)	Selection between different donor splice sites at the 5' end of an intron.	Alters the N-terminal boundary of the upstream exon.
Alternative 3' Splice Site (Alt 3'SS)	Selection between different acceptor splice sites at the 3' end of an intron.	Alters the C-terminal boundary of the downstream exon.
Mutually Exclusive Exons	Two or more exons are used in a mutually exclusive manner.	Provides structurally distinct protein variants.
Intron Retention	An intron remains in the mature mRNA instead of being spliced out.	Often leads to nonsense-mediated decay (NMD) or truncated proteins.

Statistical Power and Sample Size Considerations

Detecting differential splicing with confidence requires careful experimental design. The following table outlines key parameters influencing statistical power, drawing from large-scale benchmarking studies [29] [31].

Table 2: Parameters Influencing Statistical Power for Differential Splicing Detection

Parameter	Recommended Practice	Rationale and Impact
Biological Replicates	Minimum of 3-4 per condition; 6-8 for robust detection of subtle changes.	Accounts for natural biological variation; essential for reliable statistical testing [31].
Sequencing Depth	Typically 30-50 million paired-end reads per sample for mammalian genomes.	Ensures sufficient coverage across splice junctions for accurate PSI quantification.
Read Length	Paired-end reads (e.g., 2x100 bp or 2x150 bp) are strongly recommended.	Provides better alignment accuracy across splice junctions and isoform resolution [18].
Effect Size (ΔPSI)	Studies should be powered to detect a minimum ΔPSI of 0.1-0.2 (10-20%).	Smaller effect sizes require substantially more replicates to distinguish from technical noise [29].
Signal-to-Noise Ratio (SNR)	Aim for high SNR through quality library prep and sequencing.	Low SNR drastically reduces the ability to detect true biological differences, especially with subtle effects [29].

Computational Tools for Splicing Quantification

A wide array of specialized software tools has been developed to identify and quantify AS events from bulk RNA-Seq alignment files. The selection of a tool depends on the specific analysis goals [2] [15].

Table 3: Selected Computational Tools for Alternative Splicing Analysis from Bulk RNA-Seq

Tool Name	Core Methodology	Key Feature	Best Suited For
MAJIQ/VOILA [15]	Quantifies Local Splicing Variations (LSVs) and PSI using a Bayesian framework.	Handles complex, unannotated events; excellent for large, heterogeneous datasets.	Discovery-based analysis of complex splicing across many samples/conditions.
SUVA [28]	Detects and quantifies five types of AS events using splice junction reads.	Calculates the ratio of alternatively spliced to constitutively spliced reads (pSAR).	Focused analysis of canonical AS event types.
rMATS	Uses a hierarchical Bayesian model for PSI estimation and differential analysis.	Robust statistical model for replicate data; well-established.	Differential splicing analysis with biological replicates.
LeafCutter	Clusters intron splicing events without relying on existing transcript annotations.	Annotation-free; identifies novel intron clusters and cryptic splicing.	De novo discovery of splicing quantitative trait loci (sQTLs).
DEXSeq	Models exon usage counts to test for differential exon usage.	Generalized linear model framework; integrates well with differential expression pipelines.	Exon-level differential usage analysis.

Detailed Experimental and Computational Protocols

Protocol I: End-to-End Splicing Analysis Workflow

This protocol outlines the complete process from sample preparation to the visualization of differential splicing results.

1. Experimental Design & Sample Preparation:

• Define conditions and ensure a minimum of three biological replicates per group [31].
• Extract high-quality total RNA. Assess RNA Integrity Number (RIN > 8 recommended).
• Use stranded, paired-end library preparation protocols (e.g., Illumina TruSeq Stranded mRNA). This preserves strand information, critical for accurate transcript and splice junction assignment [18] [29].

2. Sequencing & Primary Data Generation:

• Sequence libraries to a minimum depth of 30 million paired-end reads (e.g., 2x100 bp) on an Illumina platform [29].
• Demultiplex samples and generate FASTQ files. Perform initial quality control with FastQC.

3. Bioinformatics Processing (Using nf-core/RNAseq Pipeline):

• Alignment: Use a splice-aware aligner (e.g., STAR) to map reads to the reference genome and transcriptome [18].
• Quantification: Employ alignment-based quantification with Salmon or RSEM to estimate transcript abundances. This step models uncertainty in read assignment to isoforms [18].
• Execution: The nf-core/RNAseq Nextflow pipeline automates these steps, including quality control (QC), alignment, and quantification, ensuring reproducibility [18].

4. Splicing Analysis (Using MAJIQ v2 as an Example):

• Build Splice Graph: Run the majiq build command using the genome alignment (BAM) files and a gene annotation file (GTF). This step constructs a splice graph for each gene, incorporating unannotated splice junctions [15].
• Quantify Splicing Variations: Run majiq quant on the built model to compute PSI (Ψ) values for all detected Local Splicing Variations (LSVs) in each sample [15].
• Differential Splicing Analysis: Run majiq deltapsi to compare conditions. The MAJIQ HET test is recommended for heterogeneous datasets, as it applies robust rank-based statistics without assuming a shared PSI per group [15].

5. Visualization & Interpretation:

• Use the companion tool VOILA to visualize significant LSVs. Inspect splice graphs, PSI distributions across samples, and read coverage to validate findings [15].
• Filter results for confidence (e.g., |ΔPSI| > 0.2 and probability P(|ΔPSI| > 0.2) > 0.95).
• Perform functional enrichment analysis (e.g., Gene Ontology, KEGG) on genes harboring significant differential splicing events.

Bulk RNA-Seq Splicing Analysis Workflow

Protocol II: Statistical Framework for Subtle Splicing Changes

Detecting subtle differential splicing—critical for clinical or subtype analyses—requires stringent quality control and statistical modeling [29].

1. Quality Control and Signal-to-Noise Assessment:

• Calculate the Principal Component Analysis (PCA)-based Signal-to-Noise Ratio (SNR) for your expression or PSI matrix [29].
• Interpretation: A low SNR indicates high technical noise obscuring biological signal. For studies aiming to detect subtle changes, samples or libraries with SNR < 12 should be investigated and potentially excluded [29].

2. Utilizing Spike-In Controls for Normalization:

• Spike-in controls (e.g., ERCC RNA or SIRVs) can be added during library preparation [29] [31].
• Use these controls to assess technical variability, dynamic range, and to aid in the normalization of splicing metrics, particularly when global transcriptomic changes are large.

3. Application of Heterogeneous Group Statistics:

• When analyzing datasets with high biological variability (e.g., patient samples) or combining public datasets, standard group-wise tests may fail.
• Apply the MAJIQ HET module, which uses non-parametric, rank-based tests (like Mann-Whitney U) on per-sample PSI estimates. This method is more robust to outliers and non-normal distributions common in heterogeneous data [15].

4. Batch Effect Correction:

• If samples were processed in different batches, apply batch correction methods (e.g., using limma::removeBatchEffect or ComBat) to the PSI matrix or read counts before differential testing [30].

Statistical Framework for Splicing Analysis

Research Reagent and Resource Toolkit

The following table lists essential reagents, controls, and software resources critical for ensuring the reliability and reproducibility of bulk RNA-Seq splicing profiling experiments.

Table 4: Essential Research Reagent Solutions for Splicing-Focused RNA-Seq

Category	Item	Function and Application in Splicing Analysis
Library Prep	Stranded mRNA-Seq Kit (e.g., Illumina TruSeq Stranded)	Preserves strand orientation of reads, which is mandatory for accurate annotation of splice junctions and antisense transcripts [18] [29].
Spike-In Controls	ERCC ExFold RNA Spike-In Mixes (Thermo Fisher)	A set of synthetic RNAs at known concentrations. Used to assess technical sensitivity, dynamic range, and accuracy of isoform quantification across samples [29].
Spike-In Controls	SIRV Spike-In Controls (Lexogen)	A set of synthetic isoform variants with known ratios. Specifically designed to benchmark the performance of isoform detection and quantification in RNA-Seq experiments [31].
Alignment & Quantification	STAR Aligner & Salmon	STAR performs fast, accurate splice-aware alignment to the genome [18]. Salmon performs alignment-based or alignment-free quantification, effectively modeling uncertainty in read assignment to isoforms [18].
Splicing Analysis Software	MAJIQ v2 & VOILA v2	A dedicated suite for quantifying Local Splicing Variations (LSVs), performing differential analysis (including for heterogeneous data), and visualizing complex splice graphs and results [15].
Reference Materials	Quartet Project Reference RNA Samples	Well-characterized reference materials from a family quartet. Used for inter-laboratory benchmarking, especially for assessing performance in detecting subtle differential expression and splicing [29].
Validation	TaqMan Assays or RT-qPCR Primers	Designed across specific splice junctions for independent, targeted validation of key differential splicing events identified in the RNA-Seq analysis [28].

A Practical Workflow for Splicing Analysis: From Raw Data to Biological Insight

The comprehensive analysis of alternative splicing (AS) from bulk RNA sequencing (RNA-Seq) data is a cornerstone of modern functional genomics, with direct implications for understanding disease mechanisms and identifying therapeutic targets. It is estimated that approximately 90% of human genes undergo alternative splicing, generating transcriptomic diversity crucial for development and cellular differentiation [32]. The accurate detection and quantification of splicing events—such as exon skipping, alternative 5'/3' splice site usage, and intron retention—are fundamentally dependent on the initial step of read alignment. Splice-aware aligners must precisely map sequencing reads that span exon-exon junctions to a reference genome, a computationally intensive task critical for all downstream analyses [33].

Within this framework, HISAT2 and STAR have emerged as two preeminent, yet algorithmically distinct, splice-aware aligners. The choice between them directly influences the sensitivity, accuracy, and reliability of subsequent differential splicing analysis, a key component of a thesis focused on AS in bulk RNA-Seq. HISAT2 employs a hierarchical graph FM-index (GFM) to efficiently map reads against a population of genomes, offering a balanced memory footprint [34]. In contrast, STAR utilizes a highly efficient seed-and-extend algorithm based on uncompressed suffix arrays, prioritizing mapping speed at the cost of greater memory (RAM) usage [35]. Recent studies highlight that while both tools perform robustly, subtle differences in their handling of spliced alignments can lead to variances in the identification of non-canonical junctions and variant calls, potentially impacting biological conclusions [36] [37] [33]. This article provides detailed application notes and protocols to guide researchers in selecting and implementing the optimal alignment strategy for alternative splicing research.

Comparative Analysis: HISAT2 vs. STAR for Splicing-Focused Research

Selecting the appropriate aligner requires a clear understanding of performance trade-offs. The following tables synthesize quantitative data from recent benchmarking studies to guide decision-making.

Table 1: Performance and Benchmarking Overview

Metric	HISAT2	STAR	Research Context & Implications
Core Algorithm	Hierarchical Graph FM-index (HGFM) [34]	Seed-and-extend with uncompressed suffix arrays [35]	HISAT2 is designed for efficient mapping against a graph representing genetic variation. STAR excels in raw speed for mapping to a linear reference.
Typical Mapping Rate	~93.7% - 99.5% [32] [36]	~90.5% - 99.5% [36] [35]	Both achieve high rates. HISAT2 may show marginally higher on-target rates in targeted sequencing (rLAS) [32].
Splice Junction Detection	Effective for canonical and non-canonical sites [34].	High precision for canonical and novel junctions; superior for fusion transcript detection [35].	STAR's two-pass mode is recommended for novel junction discovery [35]. Both may generate spurious junctions in repetitive regions [37].
Computational Profile	Lower memory footprint. Faster index building.	High memory usage (e.g., ~32GB for human genome). Very fast alignment [33] [35].	HISAT2 is suitable for resource-constrained environments. STAR requires high-performance hardware but reduces wall-clock time.
Impact on Downstream Variant Calling	Contributes to alignment-specific variant sets [33].	Different splice junction mapping can lead to divergent variant identification [33].	A study found less than 2% overlap in potential RNA editing sites called from different aligners, highlighting a major source of discrepancy [33].

Table 2: Decision Framework for Aligner Selection

Research Priority	Recommended Aligner	Rationale and Protocol Considerations
Standard Bulk RNA-Seq (Splicing & Expression)	STAR (if resources allow)	Superior speed for large datasets facilitates rapid iteration. Use --twopassMode Basic for comprehensive novel junction discovery [35].
Constrained Computational Resources	HISAT2	Efficient memory use allows operation on standard workstations or clusters with smaller nodes [34] [33].
Targeted or Long-Amplicon Sequencing (e.g., rLAS)	HISAT2 or STAR	Performance is comparable [32]. Choice may depend on lab pipeline familiarity. HISAT2 showed a slightly higher on-target rate in one evaluation [32].
Dual RNA-Seq (Host-Pathogen)	Combination Approach	For complex eukaryotic hosts, a pathogen-first mapping strategy is critical. Map to pathogen genome first with BWA, then unmapped reads to host genome with HISAT2/STAR to minimize misalignment [38].
Maximizing Reproducibility & Minimizing False Junctions	Either, plus EASTR post-processing	Both aligners can create false spliced alignments in repetitive regions [37]. Running the EASTR tool to filter alignments post-alignment significantly improves accuracy [37].

Detailed Experimental Protocols

The following protocols are optimized for alternative splicing analysis in a bulk RNA-Seq context.

Protocol 1: Standard Spliced Alignment Workflow with STAR for Novel Junction Discovery This protocol uses STAR's two-pass mode to maximize sensitivity for novel alternative splicing events, which is essential for discovery-based research [35].

• Genome Index Generation (One-time):
  Critical Parameters: --sjdbOverhang should be set to (Read Length - 1). Using a comprehensive annotation (GTF) during indexing improves initial junction detection [35].

• Two-Pass Alignment for Novel Splice Junction Detection:
  Application Note: The SJ.out.tab files from the first pass are pooled via --sjdbFileChrStartEnd to create a enhanced junction database for the final alignment, dramatically improving sensitivity to sample-specific splicing [35].

Protocol 2: Memory-Efficient Alignment with HISAT2 for Targeted Splicing Analysis This protocol is optimized for targeted sequencing data (like rLAS) or environments with limited RAM, where HISAT2's efficiency is advantageous [32] [34].

• Build Index with Known Splice Sites:
  Note: Building the index is separate from preparing splice site files. The index is built from the FASTA file alone [34].

• Alignment with Annotation Guidance:
  Application Note: The --dta (downstream transcriptome assembly) option reports alignments best suited for assemblers like StringTie, which is a common step in splicing isoform reconstruction. The --no-softclip option can be beneficial for targeted amplicon data where full-read alignment is desired [32] [34].

Protocol 3: Post-Alignment Filtering with EASTR to Eliminate Systematic Errors Both HISAT2 and STAR can produce false spliced alignments in repetitive regions [37]. This post-processing step is highly recommended to improve the fidelity of downstream splicing analysis.

• Run EASTR on alignment (BAM) files:

• Use filtered BAM for downstream analysis: The alignments_filtered.bam file will have alignments supporting spurious junctions removed. Studies show EASTR filtering substantially reduces false positive introns and exons during transcript assembly, improving the accuracy of novel isoform detection [37].

The Scientist's Toolkit: Essential Research Reagent Solutions

Table 3: Key Reagents, Software, and Materials for Splicing-Aware Alignment Workflows

Item Name	Type	Function in Splicing Analysis	Example/Note
rLAS (Targeted RNA Long-Amplicon Seq) Reagents [32]	Wet-lab Method	Enables deep sequencing of specific low-expressing transcripts for splicing variants.	Critical for validating splicing events in candidate disease genes cost-effectively [32].
Stranded RNA-Seq Library Prep Kit	Wet-lab Consumable	Preserves strand-of-origin information, crucial for accurately assigning reads to the correct transcript and splicing isoform.	Used in standard bulk RNA-Seq protocols for splicing.
Dual RNA-Seq Library Prep	Wet-lab Method	Allows concurrent analysis of host and pathogen transcriptomes, relevant for infectious disease splicing studies [38].	Protocol involves careful RNA isolation to preserve both host and pathogen RNA [38].
Splice-Aware Aligner (STAR/HISAT2)	Core Software	Maps RNA-seq reads across splice junctions to the reference genome.	Fundamental first step. Choice affects all downstream results [36] [33].
EASTR (Emending Alignments Tool) [37]	Post-processing Software	Identifies and removes falsely spliced alignments caused by repetitive sequences, reducing false positive splicing events.	Should be run after alignment and before transcript assembly/quantification [37].
Splicing Analysis Suite (rMATS/MAJIQ) [32]	Downstream Software	Quantifies differential splicing from aligned reads.	Performance varies; rMATS excels at exon skipping, MAJIQ handles diverse event types [32].
Reference Genome & Annotation (GTF)	Data File	Provides the coordinate system and known gene models for alignment and quantification.	Critical for alignment (--sjdbGTFfile in STAR) and quantification. Use the same version consistently.
High-Performance Compute Node	Hardware	Provides the necessary memory and CPU for alignment, especially for STAR with large genomes.	≥32 GB RAM recommended for human genome alignment with STAR [35].

Visualization of Workflows and Decision Pathways

(Diagram 1: Standard Bulk RNA-Seq Splicing Analysis Workflow)

(Diagram 2: Decision Tree for Selecting a Splice-Aware Aligner)

Event-Based vs. Isoform-Based Quantification Approaches

In bulk RNA sequencing research, the analysis of alternative splicing (AS) is fundamental for understanding transcriptome diversity, cellular differentiation, and disease mechanisms. Accurate quantification of splicing variations is computationally challenging due to the sequence similarity between isoforms originating from the same gene locus [39]. Two predominant computational paradigms have emerged to address this challenge: event-based and isoform-based quantification. Event-based methods quantify specific, localized splicing variations—such as exon skipping or intron retention—by measuring the proportion of reads supporting one event outcome over another [40]. In contrast, isoform-based methods aim to estimate the absolute or relative abundance of full-length transcript isoforms, treating each complete mRNA molecule as a quantitative unit [41].

The choice between these approaches has significant implications for a research project. Event-based analysis is often more statistically powerful for detecting specific, known changes in splicing patterns, while isoform-based analysis provides a comprehensive, systems-level view of transcriptional output but faces greater inherent uncertainty [42]. This article provides detailed application notes and protocols for both methodologies, framed within the context of a broader thesis on alternative splicing analysis using bulk RNA-seq. It is designed to equip researchers, scientists, and drug development professionals with the knowledge to select, implement, and interpret these critical analytical strategies.

Event-based and isoform-based quantification represent fundamentally different ways of modeling the transcriptome. Their core differences in methodology directly influence their applications, strengths, and limitations.

Event-Based Quantification focuses on discrete, categorical splicing changes. It reduces the complexity of the transcriptome by examining predefined types of AS events, such as Skipped Exons (SE), Alternative 5' or 3' Splice Sites (A5SS, A3SS), Mutually Exclusive Exons (MXE), and Retained Introns (RI) [40]. The output is typically a "Percent Spliced In" (PSI or Ψ) value, representing the proportion of transcripts from a gene that include a particular sequence element (e.g., an exon). This approach is highly interpretable, as it links directly to specific mechanistic changes in RNA processing. Tools like rMATS, MAJIQ, and SUPPA2 are prominent examples [42] [43].

Isoform-Based Quantification attempts to solve the more general problem of estimating the abundance of every annotated (and sometimes novel) transcript isoform. Methods like Salmon, Kallisto, and RSEM use probabilistic models to resolve the origin of sequencing reads that map to multiple isoforms, outputting estimates in counts or Transcripts Per Million (TPM) [39] [41]. This provides a complete picture of transcriptional output but is notoriously difficult due to the shared exonic sequences among isoforms, especially with short-read data.

The table below summarizes the key characteristics of the two approaches:

Table 1: Core Methodological Comparison of Quantification Approaches

Aspect	Event-Based Quantification	Isoform-Based Quantification
Quantitative Unit	Local splicing event (e.g., exon, splice site).	Full-length transcript isoform.
Typical Output Metric	Percent Spliced In (PSI/Ψ).	Counts, TPM, or FPKM.
Primary Challenge	Accurate read counting for event junctions and bodies.	Resolving multi-mapping reads to correct isoforms.
Interpretability	Directly links to a specific biological mechanism.	Requires downstream analysis to infer splicing changes.
Examples of Tools	rMATS, MAJIQ, SUPPA2, DEXSeq [42] [40].	Salmon, Kallisto, RSEM, StringTie [39] [41].
Best Suited For	Detecting specific, known types of differential splicing.	Global profiling, isoform switching, novel isoform discovery.

Performance benchmarks indicate that the choice of approach and tool significantly impacts results. In a systematic evaluation of differential splicing tools, exon-based methods (a subtype of event-based) and certain event-based tools like rMATS and MAJIQ generally scored well on precision and recall [42] [43]. Meanwhile, studies on isoform quantification show that while alignment-free tools like Salmon and Kallisto are fast and accurate on idealized data, accuracy degrades with complex gene structures, lowly expressed transcripts, and incomplete annotations [39] [41]. Notably, isoform-based quantification of genes like TP53, which produces isoforms with antagonistic functions (e.g., pro-apoptotic full-length vs. anti-apoptotic Δ133 variants), is critical for understanding functional outcomes in cancer [39].

Diagram 1: Conceptual workflow for event and isoform analysis.

Detailed Experimental Protocols

Protocol for Event-Based Differential Splicing Analysis with rMATS

This protocol uses rMATS (replicate Multivariate Analysis of Transcript Splicing), a widely used, computationally efficient tool for detecting differential splicing from replicate RNA-seq experiments [40].

1. Prerequisite Data Preparation:

• Sequence Alignment: Align RNA-seq reads (paired-end recommended) to the reference genome using a splice-aware aligner such as STAR [18]. Output coordinate-sorted BAM files for each sample.
• Sample Grouping: Create a text file listing the BAM file paths, separated into two groups (e.g., treatment vs. control).

2. Running rMATS Quantification: Execute rMATS using the following command structure:

Critical Parameters:

• --b1/--b2: Text files containing paths to BAM files for each condition.
• --gtf: Reference transcriptome annotation in GTF format.
• --readLength: Must match your sequencing read length.
• --libType: Strandedness of your library (fr-unstranded, fr-firststrand, fr-secondstrand).
• --novelSS: Enable to detect splicing events using novel (unannotated) splice sites.

3. Output Interpretation: rMATS generates results for five event types: Skipped Exon (SE), Alternative 5' Splice Site (A5SS), Alternative 3' Splice Site (A3SS), Mutually Exclusive Exon (MXE), and Retained Intron (RI). The key output file JC.raw.input.*.txt contains:

• Inclusion/Exclusion Counts: Read counts supporting the included and excluded isoforms.
• Inclusion Level (Ψ): Calculated for each sample group.
• P-value and FDR: Statistical significance of the difference in Ψ between groups.
• A significant splicing change is typically defined by |ΔΨ| > a threshold (e.g., 0.1 or 0.2) and FDR < 0.05.

Protocol for Isoform-Based Quantification and Analysis with Salmon

This protocol employs Salmon, a rapid, alignment-free tool that performs transcript-level quantification via quasi-mapping [39]. It is part of a best-practice workflow that also facilitates quality control [18].

1. Indexing the Transcriptome: First, build a Salmon index from a reference transcriptome (FASTA file derived from a genome and GTF annotation).

2. Quantification of Samples: Run Salmon in quasi-mapping mode on each sample's raw FASTQ files.

Critical Parameters:

• -l A: Lets Salmon automatically infer the library type.
• --gcBias: Corrects for GC-content bias.
• --validateMappings: Enables a more accurate, validation-based mapping strategy.
• For maximum quality control integration, a hybrid workflow using STAR for alignment followed by Salmon in alignment-based mode (-a flag) is recommended [18].

3. Aggregating Results for Differential Analysis: Use the tximport package in R/Bioconductor to import Salmon's transcript-level abundance estimates (quant.sf files) and aggregate them to the gene level, correctly handling abundance estimation uncertainty.

The resulting object txi contains gene-level count estimates suitable for differential expression analysis with tools like DESeq2 or limma-voom. For differential transcript usage (DTU) analysis, packages like DRIMSeq or DEXSeq (applied at the transcript level) can be used on the transcript-level counts.

The Scientist’s Toolkit: Essential Reagents and Materials

Table 2: Key Research Reagent Solutions for Splicing Analysis

Reagent/Material	Function in Analysis	Considerations for Selection
Reference Genome & Annotation (e.g., GENCODE, Ensembl)	Provides the coordinate system and transcript models for read alignment and quantification.	Use the most comprehensive, version-matched annotation available (e.g., GENCODE). Incompleteness severely impacts isoform-based methods [41].
Splice-Aware Aligner (e.g., STAR)	Aligns RNA-seq reads across exon junctions to the genome, producing BAM files essential for most event-based tools and QC [18].	Chosen for speed and accuracy. Alignment is a prerequisite for rMATS and optional but recommended for Salmon-based workflows.
Quantification Software Suite	Core analytical engine (e.g., rMATS, Salmon, Kallisto).	Choice dictates the paradigm. For event-based analysis, rMATS is computationally efficient [40]. For isoform-based, alignment-free tools like Salmon are fast and accurate [39].
High-Quality, Strand-Specific RNA-seq Library	The primary source of data. Library preparation protocol influences analysis.	Paired-end, strand-specific libraries are strongly recommended [18]. They provide more information for resolving transcript origins and correct strand assignment of splicing events.
Synthetic RNA Spike-Ins (e.g., RNA Sequins)	Act as internal controls for benchmarking isoform detection and quantification accuracy in complex samples [44].	Used to assess the technical performance of the wet-lab and computational pipeline, especially for novel isoform detection or method benchmarking.
Differential Expression/Usage Analysis Package (e.g., DESeq2, limma, DRIMSeq)	Performs statistical testing for differential gene expression (DGE) or differential transcript usage (DTU) from count matrices.	DESeq2/limma are standard for DGE. For DTU from isoform counts, specialized tools like DRIMSeq or the spliceDT method in DEXSeq are required.

Performance Characteristics and Practical Considerations

The practical application of these methods requires an understanding of their performance under different experimental conditions. Systematic evaluations provide critical guidance.

Table 3: Performance and Computational Characteristics of Representative Tools

Tool	Category	Key Performance Insight	Computational Demand	Best Use Case
rMATS [42] [40]	Event-Based	High precision for SE, A5SS, A3SS events; lower concordance on RI events.	Moderate job time, efficient CPU/RAM use. Replicate number exponentially increases time [40].	Differential splicing analysis with biological replicates.
Salmon [39] [45]	Isoform-Based (Alignment-free)	High accuracy (correlation >0.9) on idealized data; accuracy decreases with transcript complexity & short length.	Very fast; minimal RAM in quasi-mapping mode.	Large-scale studies, rapid isoform quantification.
Kallisto [39] [45]	Isoform-Based (Alignment-free)	Performance comparable to Salmon. High accuracy in simulations for full-length scRNA-seq protocols [45].	Very fast; low memory footprint.	Similar to Salmon; often used in single-cell workflows.
MAJIQ [42]	Event-Based	Scored well on precision/recall benchmarks. Provides local splicing quantification.	Moderate to high.	Analysis where quantifying PSI at single-nucleotide resolution is needed.
RSEM [39] [41]	Isoform-Based (Alignment-based)	Good accuracy, but often slower. Performance stable across simulation types.	Higher than alignment-free tools; requires pre-aligned BAMs.	When integrated into workflows relying on genome-aligned BAMs for QC.

Critical Decision Factors:

• Biological Question: For focused hypothesis testing on specific splicing mechanisms (e.g., exon skipping in cancer), event-based methods are direct and powerful. For exploratory discovery of isoform switching—where the dominant transcript from a gene changes between conditions—isoform-based quantification followed by DTU analysis is necessary [39].
• Annotation Completeness: Isoform-based methods are heavily dependent on a complete transcript annotation. Performance degrades if the expressed isoforms are not in the reference [41]. Event-based methods can more easily accommodate unannotated splice sites (e.g., with rMATS --novelSS flag).
• Sequencing Technology: Standard short-read RNA-seq necessitates the indirect methods described. The advent of long-read sequencing (PacBio, Oxford Nanopore) mitigates the multi-mapping problem, allowing for direct observation of full-length isoforms. Tools like IsoQuant, Bambu, and StringTie2 are benchmarked for long-read isoform detection [44]. While not yet standard for quantification due to cost and error rates, long-read data can be used to refine annotations for subsequent short-read analysis.
• Single-Cell RNA-seq: Most high-throughput scRNA-seq protocols (e.g., 10x Genomics) have strong 3'-end bias, making full-length isoform quantification extremely challenging. Specialized tools like Scasa [46] and SCALPEL [47] have been developed to address this, often showing that isoform-level analysis can reveal cell subtypes invisible to gene-level analysis. However, bulk analysis remains the primary and more reliable method for splicing studies.

Diagram 2: A unified analytical pipeline for RNA-seq splicing analysis.

Event-based and isoform-based quantification are complementary frameworks for alternative splicing analysis in bulk RNA-seq. The optimal choice is not universal but depends on the specific research objectives, available genomic resources, and experimental design.

For a typical thesis project in bulk RNA-seq research, the following strategic approach is recommended:

• Define the Primary Objective: If the goal is to identify specific changes in known splicing patterns (e.g., in a disease model), begin with an event-based analysis using rMATS. It is computationally efficient, statistically robust for replicates, and yields easily interpretable results [42] [40].
• Employ a Hybrid Workflow for Discovery: For a comprehensive transcriptomic profile, implement a Salmon-based isoform quantification pipeline as the primary analysis [18]. Use tximport and DESeq2 to perform differential gene expression and DRIMSeq for differential transcript usage. This can reveal global isoform switching events.
• Validate and Integrate Findings: Significant events from both analyses should be validated (e.g., by RT-PCR). Importantly, cross-reference the results: a significant exon skip (from rMATS) should correspond to a change in the abundance of isoforms differing by that exon (from Salmon).
• Acknowledge Technological Limitations: Recognize that short-read RNA-seq has inherent limitations in resolving highly similar isoforms. For critical genes or when novel isoforms are suspected, consider targeted validation with long-read sequencing or orthogonal assays [44].

Ultimately, the convergence of evidence from both event-centric and isoform-centric viewpoints will provide the most robust and biologically insightful conclusions for any thesis on alternative splicing.

Alternative splicing (AS) is a fundamental post-transcriptional process enabling tremendous transcriptomic and proteomic diversity from a limited set of genes. In humans, over 90-95% of multi-exonic genes undergo alternative splicing, making its study critical for understanding cellular function, development, and disease mechanisms such as cancer and neurodegeneration [48] [49]. Bulk RNA sequencing (RNA-seq) has become the standard technology for genome-wide profiling of splicing variations, but the short nature of sequencing reads (typically 100-150 bp) poses challenges for accurately discerning between highly similar isoforms [48].

Computational tools for differential splicing analysis have evolved to address these challenges, employing different methodological frameworks. These can be broadly categorized into exon-based, event-based, and isoform-based methods [48]. This article focuses on four leading tools—MAJIQ, rMATS, LeafCutter, and DEXSeq—detailing their applications, protocols, and integration within a comprehensive thesis research framework for bulk RNA-seq analysis. The selection of an appropriate tool is not trivial, as a systematic evaluation has shown that performance varies considerably across different datasets and sample sizes [48].

Table 1: Methodological Overview and Classification of Featured Tools

Tool	Core Method	Primary Unit of Analysis	Key Splicing Events Detected	Statistical Foundation
DEXSeq [48]	Exon-based	Individual exonic bins	Differential exon usage	Generalized linear model (Negative Binomial) [49]
rMATS [50]	Event-based	Pre-defined canonical events (SE, RI, etc.)	Exon skipping, intron retention, alternative splice sites	Hierarchical Bayesian model [50]
LeafCutter [51]	Intron-based	Clusters of alternative introns	Any variation altering intron excision	Dirichlet-Multinomial generalized linear model [52]
MAJIQ [15]	Event-based (LSVs)	Local Splicing Variations (LSVs)	Complex and de novo junctions, including unannotated events	Bayesian inference with non-parametric tests [15]

DEXSeq: Detailed Exon Analysis

DEXSeq (Differential Exon Usage analysis) is an exon-based method that operates by testing for deviations in read counts of individual exonic regions relative to the overall gene expression [48] [49]. It extends the statistical framework of DESeq2, using a negative binomial generalized linear model (GLM) to account for biological variability and over-dispersion in read counts [49]. Within a thesis on alternative splicing, DEXSeq is particularly powerful for hypothesis-driven research where the focus is on identifying specific exons that show differential usage between conditions, such as in studies of tissue-specific regulation or the functional impact of mutations in exonic enhancer/silencer elements.

Experimental Protocol

The typical workflow for DEXSeq analysis, as implemented in R/Bioconductor, involves the following key steps [48] [49]:

• Input Preparation: Generate count data for non-overlapping exonic bins. This is often done using the DEXSeq Python script (dexseq_count.py) which takes aligned reads in BAM/SAM format and a gene annotation file in GFF format. The annotation is flattened into a set of exonic counting bins.

• Data Import and Normalization in R: Read the count files into an ExonCountSet object, which holds counts, gene IDs, exon bin IDs, and sample metadata.

• Statistical Testing and Result Extraction: Fit the GLM and test for differential exon usage. The core test assesses whether the interaction term condition:exon is significant, indicating that the exon's usage differs between conditions relative to the gene total.

• Visualization: Generate plots for interpretation, such as the per-exon expression or splicing heatmaps.

Output Interpretation

The primary result object (DEXSeqResults) contains rows for each tested exonic bin. Key columns include:

• groupID: The gene identifier.
• featureID: The exon/bin identifier.
• pvalue and padj: The raw and adjusted p-values for the significance of differential usage.
• log2fold: The log2 fold change in the exon's usage between conditions. A significant adjusted p-value (e.g., padj < 0.05) indicates that the proportion of reads mapping to that specific exonic bin has changed significantly relative to the overall expression of the gene [49].

Advantages and Limitations

Advantages: DEXSeq is highly sensitive for detecting differential exon usage, benefits from the robust statistical framework of DESeq2, and is well-integrated into the Bioconductor ecosystem for downstream analysis [48] [49]. Limitations: Its focus on exonic bins can make biological interpretation of specific splicing events (like exon skipping) less direct compared to event-based methods. It also relies heavily on the accuracy of the provided gene annotation [48] [53].

DEXSeq Analysis Workflow: From aligned reads to significant exons.

rMATS: Robust Analysis of Canonical Splicing Events

rMATS (replicate Multivariate Analysis of Transcript Splicing) is a widely used event-based tool designed for the robust detection of differential alternative splicing from RNA-seq data with biological replicates [50]. Its statistical strength lies in a hierarchical model that simultaneously accounts for estimation uncertainty within individual replicates (due to sequencing depth) and biological variability among replicates [50]. For thesis research, rMATS is exceptionally valuable for case-control studies aiming to quantify changes in well-defined, canonical splicing event types, such as identifying exon skipping events in cancer versus normal tissues or intron retention under cellular stress. It also uniquely supports the analysis of paired replicate designs (e.g., matched tumor-normal pairs), enhancing statistical power in clinical study designs [50] [54].

Experimental Protocol

Running rMATS typically involves a single command that handles alignment (if starting from FASTQ) or quantification (if starting from BAM) [54].

• Input Preparation: Prepare two text files (s1.txt, s2.txt) listing the paths to replicate files for each condition. For FASTQ inputs, the format within the file is rep1_R1.fq:rep1_R2.fq,rep2_R1.fq:rep2_R2.fq. BAM files can also be used directly.

• Execution: The main analysis is run via the rmats.py script. Key parameters must be specified.

  • --s1/--s2: Sample lists for groups 1 and 2.
  • --gtf: Reference transcriptome annotation.
  • --bi: Path to a STAR genome index (required for FASTQ input).
  • --od: Output directory.
  • --cstat: Cutoff for splicing difference (∆PSI). Default 0.0001 (0.01%) [54].
  • --tstat: Number of threads for statistical modeling.

• Output Interpretation: rMATS produces separate output files for each event type: Skipped Exon (SE), Retained Intron (RI), Alternative 5' Splice Site (A5SS), Alternative 3' Splice Site (A3SS), and Mutually Exclusive Exons (MXE). The most commonly used files are the .MATS.JC.txt (junction counts only) or .MATS.JCEC.txt (junction counts + exon body counts). Critical columns include:

  • IncLevel1, IncLevel2: Mean inclusion level (PSI or Ψ) for each group.
  • IncLevelDifference: ∆PSI (PSIgroup2 - PSIgroup1).
  • PValue, FDR: Significance of the differential splicing event. A significant event is typically defined by FDR < 0.05 and |IncLevelDifference| > a threshold (e.g., 0.05 or 0.1) [50].

Advantages and Limitations

Advantages: rMATS is highly robust for replicate data, offers flexible hypothesis testing on the magnitude of ∆PSI, and is the standard tool for analyzing canonical splicing events. Its paired replicate model is a unique strength for clinical studies [48] [50]. Limitations: It is limited to five pre-defined event types and may miss complex or unannotated splicing variations. The analysis can be computationally intensive for large datasets [15] [55].

rMATS Splicing Event Detection and Classification.

LeafCutter: Annotation-Free Intron-Centric Analysis

LeafCutter adopts a novel, annotation-free approach that quantifies splicing variation through differential intron excision rather than exon inclusion [51]. It clusters together introns that share a splice site, defining "intron clusters" where the relative usage of each intron is quantified. This allows it to detect any splicing variation that alters an intron boundary, including complex and unannotated events, without relying on a pre-defined transcriptome [51] [55]. In thesis research, LeafCutter is particularly powerful for exploratory analysis in non-model organisms, for identifying novel splicing events in disease (like cancer), or for large-scale studies like splicing quantitative trait locus (sQTL) mapping due to its computational efficiency [51] [52]. Its variant, LeafCutterMD, is specifically designed for detecting outlier splicing events in individual samples, making it highly relevant for rare Mendelian disease research [52].

Experimental Protocol

The LeafCutter workflow involves separate steps for preparation and statistical analysis.

• Intron Quantification: Generate junction read counts from BAM files using the leafcutter_cluster_regtools.py script (or the regtools suite).

  This generates two key files: *_perind_numers.counts.gz (junction counts) and *_perind.counts.gz (intron cluster definitions).

• Differential Splicing Analysis in R: Use the LeafCutter R package to perform the group comparison.

• Visualization: LeafCutter provides functions to visualize the differential intron clusters.

Output Interpretation

The differential splicing results highlight intron clusters that show significant changes in the relative usage of their constituent introns. The primary metric is the intron excision ratio. A significant result indicates a shift in which introns within a cluster are being excised (spliced out), pointing to a change in splicing patterns for that genomic region. For LeafCutterMD, the output identifies specific sample-intron pairs where the intron usage is a statistical outlier compared to the cohort, pinpointing potential aberrant splicing in individual patients [52].

Advantages and Limitations

Advantages: Annotation-free, enabling discovery of novel events; captures complex splicing variations; computationally efficient and scalable to thousands of samples [51] [55]. Limitations: Biological interpretation can be challenging, as results are in terms of intron clusters rather than traditional event types. It requires sufficient junction coverage for reliable quantification [15] [53].

MAJIQ: Modeling Complex and De Novo Splicing Variations

MAJIQ (Modeling Alternative Junction Inclusion Quantification) is an advanced event-based framework that introduces the concept of Local Splicing Variations (LSVs) [15]. An LSV is defined as any splice graph node (an exon or intron) with more than one incoming or outgoing junction. This formulation allows MAJIQ to capture not only canonical splicing events but also complex variations involving multiple junctions and de novo (unannotated) splice sites [15] [56]. For a thesis involving deep investigation of splicing regulation, MAJIQ is invaluable for analyzing large, heterogeneous datasets (e.g., from consortia like GTEx), studying tissue-specific splicing complexity, or identifying disease-specific aberrant splicing that deviates from known annotations [15]. Its accompanying VOILA visualization tool provides unparalleled interactive exploration of splice graphs.

Experimental Protocol

The MAJIQ v2 pipeline consists of distinct modules for building, quantifying, and visualizing [15] [56].

• Builder: Constructs the splice graph from BAM files and annotation (GFF3), incorporating de novo junctions.

  The config.txt file lists BAM files and their groups.

• Quantifier: Quantifies LSVs and computes differential splicing (dPSI) between conditions. MAJIQ employs a Bayesian inference model to estimate Posterior Distributions (P(Ψ) or P(ΔΨ)).

• Visualization with VOILA: Generates interactive HTML files for exploring LSVs, their quantification, and splice graphs.

Output Interpretation

MAJIQ's core output includes the probability that the absolute change in splicing (|ΔPSI|) exceeds a user-defined threshold (e.g., P(|ΔPSI| > 0.1 or 0.2)) [15]. Analysts can filter LSVs based on this confidence measure and the expected ΔPSI value. The VOILA viewer allows direct inspection of the splice graph, showing the read coverage, junction counts, and PSI values for each sample group, which is critical for interpreting complex LSVs.

Advantages and Limitations

Advantages: Uniquely captures complex and unannotated splicing events; powerful Bayesian quantification with confidence intervals; exceptional visualization for result interpretation; designed for large datasets [15] [56]. Limitations: Has a steeper learning curve due to its conceptual complexity. The analysis can be memory-intensive for very large genomes or deep sequencing data [55].

MAJIQ v2 Pipeline: From data building to interactive visualization.

Table 2: Research Reagent Solutions for Differential Splicing Analysis

Reagent / Resource	Function in Analysis	Example/Tool	Critical Notes
Reference Genome & Annotation	Provides genomic coordinates for alignment and feature identification.	ENSEMBL GTF/GFF3 files, GENCODE	Version consistency between alignment, quantification, and annotation is crucial [54].
Alignment Software	Maps RNA-seq reads to the reference genome, identifying splice junctions.	STAR, HISAT2	STAR is preferred for sensitivity in junction detection [55] [54].
Quantification Scripts	Generate count data for specific features (exons, junctions, introns).	dexseq_count.py (DEXSeq), regtools (LeafCutter)	Must match the counting strategy of the downstream differential tool [49].
Statistical Software Environment	Provides the computational framework for running analyses.	R/Bioconductor (DEXSeq), Python (rMATS, MAJIQ, LeafCutter)	Manage dependencies via Conda or Docker for reproducibility [54].
High-Performance Computing (HPC) Resources	Provides necessary CPU, memory, and storage for data-intensive steps.	Cluster with SLURM/SGE scheduler	rMATS and whole-genome alignment are particularly resource-intensive [54].

Selecting the optimal tool from MAJIQ, rMATS, LeafCutter, and DEXSeq is a consequential decision that shapes the scope and findings of thesis research. The choice should be driven by the specific biological question, the nature of the RNA-seq data, and the desired balance between discovery and precision.

For a thesis focused on novel splicing discovery in heterogeneous samples or large cohorts (e.g., across diverse tissues or patient populations), MAJIQ and LeafCutter are superior due to their ability to capture unannotated and complex events [15] [51]. If the research aims to precisely quantify changes in canonical splicing events (like exon skipping) in a controlled case-control study with replicates, rMATS offers a robust, statistically rigorous framework [48] [50]. When the hypothesis targets differential usage of known exonic regions, perhaps in conjunction with differential gene expression analysis, DEXSeq provides a sensitive and integrated approach within the Bioconductor ecosystem [49] [53].

A comprehensive thesis may strategically employ multiple tools: using LeafCutter for an unbiased, genome-wide survey of splicing changes, followed by rMATS for precise quantification of high-interest canonical events, and finally MAJIQ's VOILA for deep visualization and validation of complex candidates. This multi-tool strategy, framed by a clear understanding of each method's strengths, will yield a more robust and insightful analysis of alternative splicing, forming a solid computational foundation for a graduate thesis in modern transcriptomics.

Core Quantitative Metrics and Computational Pipelines

The analysis of alternative splicing (AS) from bulk RNA sequencing data centers on the accurate quantification of splicing changes between biological conditions. The principal metric for this is Percent Spliced In (PSI or Ψ), which quantifies the relative inclusion of a specific splice junction or exon. For a given splicing event, PSI is calculated as the ratio of reads supporting the inclusion of a feature to the total reads covering that event [15]. The differential change in PSI between conditions, denoted as dPSI (ΔΨ), is the key measure of biologically significant splicing alteration, typically considered meaningful when |ΔΨ| > 0.1 (or 10%) [15].

Large-scale and heterogeneous RNA-seq datasets pose significant challenges for splicing quantification, including increased technical variability and the prevalence of unannotated splice variants [15]. Modern computational pipelines like MAJIQ v2 are designed to address these issues. This pipeline employs a Bayesian framework to model the posterior distribution of PSI and ΔΨ for complex Local Splicing Variations (LSVs), accounting for read distribution and coverage [15]. For heterogeneous data, its MAJIQ HET module uses robust, non-parametric rank-based statistics (e.g., Mann-Whitney U test) on per-sample PSI estimates, improving power and reproducibility [15].

Table 1: Key Bioinformatics Tools for PSI/dPSI Quantification from Bulk RNA-Seq

Tool Name	Statistical Core	Quantification Unit	Key Strength	Best Suited For
MAJIQ v2 [15]	Bayesian, Non-parametric (HET)	Local Splicing Variations (LSVs)	Handles complex/unannotated events, large heterogeneous datasets	Discovery analysis in large cohorts (e.g., GTEx, TCGA)
rMATS [49]	Likelihood-ratio test	Pre-defined splicing event types (SE, RI, etc.)	High sensitivity for canonical events, widely used	Controlled experiments with replicates, focused hypothesis testing
DEXSeq [49]	Generalized Linear Model (Negative Binomial)	Exonic bins	Tests for differential exon usage, integrates with DESeq2 workflow	Condition-specific exon inclusion, paired with gene expression analysis
LeafCutter [49]	Non-parametric	Intron excision clusters	Event-free, identifies novel splice junctions from intron boundaries	Unbiased discovery of splicing quantitative trait loci (sQTLs)

The choice of tool dictates the analytical workflow. A standard bulk RNA-seq splicing analysis begins with read alignment (using STAR or HISAT2) and the generation of splice junction maps [49]. The selected quantification tool (e.g., MAJIQ builder) then constructs a splice graph, identifies splicing events, and calculates PSI per sample [15]. Finally, statistical comparison between defined experimental groups yields dPSI values and confidence metrics (e.g., P(|ΔΨ| > 0.1)) for downstream biological interpretation [15].

Experimental Protocols for Splicing Analysis

Protocol 1: Bulk RNA-Seq for Differential Splicing Analysis

This protocol outlines the steps from tissue to quantified splicing events, optimized for detection of dPSI.

• Sample Preparation & RNA Isolation: Homogenize tissue or cell pellets in TRIzol. Isolate total RNA following standard phase-separation protocols. Assess RNA integrity using an Agilent Bioanalyzer (RIN > 8.0 required).
• Library Preparation: Deplete ribosomal RNA using species-specific probes (e.g., Illumina Ribo-Zero Plus). Prepare stranded RNA-seq libraries using a kit such as Illumina TruSeq Stranded Total RNA. Use poly-A selection only if focusing on mature mRNA.
• Sequencing: Perform paired-end sequencing (2x150 bp) on an Illumina NovaSeq platform. Target a minimum depth of 50-100 million reads per sample to ensure sufficient coverage of splice junctions [57].
• Bioinformatic Processing:
  • Quality Control & Trimming: Use FastQC for quality assessment and Trimmomatic to remove adapters and low-quality bases.
  • Alignment: Align reads to the reference genome (e.g., GRCh38) using the splice-aware aligner STAR [49]. Generate alignment files in BAM format, sorted and indexed.
  • Splice Junction Quantification: Process BAM files with the MAJIQ builder to construct splice graphs and quantify LSVs [15]. Alternatively, generate junction count files for input into rMATS or LeafCutter.
• Differential Splicing Analysis: Run the MAJIQ quantifier to calculate PSI per sample and dPSI between experimental groups. Apply a filter of |ΔΨ| > 0.1 and a probability P(|ΔΨ| > 0.1) > 0.95 for high-confidence events [15]. Visualize results using VOILA or custom sashimi plots.

Protocol 2: Integrated Long-Read and Short-Read Validation

Long-read sequencing (PacBio or Nanopore) captures full-length isoforms, providing orthogonal validation for short-read-based dPSI calls and discovering novel isoforms [57].

• Parallel Library Preparation: From the same RNA aliquot used for short-read sequencing, prepare an independent library for long-read sequencing.
  • For PacBio Iso-Seq: Perform reverse transcription with a template-switching oligo to create full-length cDNA. Use size selection to enrich for transcripts >2 kb. Prepare SMRTbell libraries and sequence on the Sequel IIe system [57].
  • For Nanopore: Use the cDNA-PCR sequencing kit (SQK-PCS112). Perform PCR with barcoding for multiplexing. Sequence on a PromethION flow cell.
• Long-Read Data Analysis:
  • PacBio: Process subreads to generate circular consensus sequences (CCS). Classify full-length non-chimeric reads and cluster them into isoforms using the Iso-Seq3 pipeline. Annotate isoforms against a reference using SQANTI3 [57].
  • Nanopore: Basecall with Guppy, align with minimap2, and collapse duplicates. Perform isoform identification with tools like FLAIR or StringTie2.
• Integration & Validation: Map the long-read-derived isoforms back to the genome. Compare the exon-intron structure of these full-length isoforms to the splicing events (LSVs) identified by MAJIQ from short-read data. Confirm that changes in the abundance of specific full-length isoforms align with the direction and magnitude of predicted dPSI values.

Case Study: Splicing Analysis in Neuronal Reprogramming Research

A critical application of dPSI quantification is in mechanistic studies of cellular reprogramming, such as the debated conversion of astrocytes to neurons upon depletion of the RNA-binding protein PTBP1 [58].

Experimental Design: Researchers used an inducible, astrocyte-specific Ptbp1 knockout mouse model (Aldh1l1-Cre/ERT2; Ptbp1^loxP/loxP) with a lineage tracer. This allowed for the precise isolation of genetically labeled astrocytes lacking PTBP1 [58].

Bulk RNA-seq & Splicing Quantification: RNA was extracted from fluorescence-activated cell-sorted (FACS) lineage-traced astrocytes from knockout and control animals. Standard bulk RNA-seq libraries were prepared and sequenced. The data was analyzed using a splicing quantification pipeline (e.g., MAJIQ or rMATS) to compute dPSI for all detected splicing events between knockout and control groups [58].

Key Quantitative Findings:

• PTBP1 depletion caused widespread splicing alterations in astrocytes, with hundreds of significant dPSI events [58].
• Crucially, the pattern of these dPSI changes was distinct from the canonical neuronal splicing pattern. The ΔΨ values for neuronal-specific exons in reprogrammed astrocytes did not match the ΔΨ values observed in genuine neurons [58].
• This quantitative dPSI analysis provided direct molecular evidence that PTBP1 loss alone does not install a neuronal splicing program, challenging the reprogramming hypothesis and highlighting the independent role of PTBP1 in regulating astrocyte-specific splicing [58].

Table 2: Summary of Splicing Analysis Workflows and Applications

Workflow Type	Core Technology	Typical Output Metrics	Key Application in Thesis Research	Considerations
Standard Bulk Analysis	Short-read Illumina RNA-seq	PSI, dPSI, P(	ΔΨ	>threshold)	Hypothesis-driven: Test effect of a perturbation (e.g., gene knockout, drug) on splicing [58].	Requires high coverage; result is population average.
Discovery in Cohorts	Heterogeneous bulk RNA-seq datasets [15]	ΔΨ distributions, novel LSV identification	Observational: Identify splicing biomarkers in disease cohorts or across tissues [15].	Needs methods robust to batch effects (e.g., MAJIQ HET).
Isoform Resolution	Long-read sequencing (PacBio/Nanopore) [57]	Full-length isoform sequences, counts	Validation & Discovery: Resolve complex loci, validate short-read dPSI, discover novel isoforms [57].	Higher cost, lower throughput; often used on subset of samples.
Spatial Context	Visium Spatial Transcriptomics (long-read) [59]	PSI per spatial spot	Spatial Regulation: Map splicing gradients or changes in tissue morphology (e.g., cortical layers) [59].	Emerging technology; lower throughput and higher computational cost.

Table 3: Key Research Reagent Solutions for Splicing Analysis

Reagent / Resource	Function	Example Product / Tool
Ribosomal RNA Depletion Kit	Removes abundant rRNA to increase sequencing depth of mRNA and non-coding RNA, crucial for splicing junction coverage.	Illumina Ribo-Zero Plus, NEBNext rRNA Depletion Kit
Stranded RNA-seq Library Prep Kit	Preserves strand information during cDNA synthesis, allowing accurate assignment of reads to the correct transcriptional strand and splice junction orientation.	Illumina TruSeq Stranded Total RNA, KAPA RNA HyperPrep
Splice-Aware Aligner	Software that accurately maps RNA-seq reads across splice junctions, essential for building correct splice graphs.	STAR [49], HISAT2 [49]
Splicing Quantification Software	Core tool for calculating PSI and dPSI from aligned RNA-seq reads. Choice depends on experimental design.	MAJIQ [15], rMATS [49], DEXSeq [49]
Long-read Sequencing Kit	Generates full-length cDNA sequences for isoform validation and discovery, complementing short-read dPSI data.	PacBio Iso-Seq HiFi Kit, Oxford Nanopore cDNA-PCR Sequencing Kit
Splicing Factor Modulators	Small molecules or ASOs used as perturbagens to experimentally alter splicing for functional validation [60].	Risdiplam (SMN2 splicing modifier), PTBP1-targeting ASOs [58]

Visualization of Analysis Workflows

Bulk RNA-Seq Splicing Analysis Pipeline [15] [49]

Integrated Short-Read and Long-Read Splicing Analysis [57]

The Quantitative Significance of Alternative Splicing in Transcriptomics

Systematic analyses of bulk RNA-seq datasets establish that alternative splicing is a major, independent regulator of cellular function, contributing nearly equally with gene expression changes to the observed biological signal. Investigations across hundreds of diverse human RNA-seq datasets reveal that splicing is not a peripheral mechanism but is central to adaptive cellular responses.

Table 1: Prevalence of Differential Splicing vs. Differential Expression [61]

Metric	Differential Expression (Average)	Differential Splicing (Average)	Overlap (Average)
Number of Significant Genes	4,327 genes (26.1% of tested)	2,247 genes (12.9% of tested)	1,252 genes
Proportion of Signal	51.8% of combined biological signal	48.1% of combined biological signal	Not applicable
Gene-Level Overlap	33.6% of DE genes are also differentially spliced	55.7% of DS genes are also differentially expressed	Shared regulatory node

The biological signals conveyed by expression and splicing are distinct yet intertwined. Gene-set enrichment analysis (GSEA) demonstrates that while some pathways are co-regulated, many are uniquely modulated by one mechanism [61]. For instance, in a study of TGFβ-treated epithelial cells, telomere maintenance pathways were enriched exclusively among differentially spliced genes, not differentially expressed genes, highlighting a unique regulatory layer [61]. The median correlation between enrichment scores for splicing and expression changes across pathways is low (Spearman’s ρ = 0.26), confirming they often capture independent biological information [61].

Table 2: Contribution of Splicing to Gene Function and Analysis Outcomes [61]

Observation	Quantitative Finding	Implication
Isoform Expression Contribution	Changing isoforms contribute 45.6% of parent gene expression on average.	Splicing changes are not artifacts of lowly expressed isoforms but affect central gene products.
Pervasiveness of Splicing	93% of multi-isoform genes are differentially spliced in at least one analyzed condition.	Splicing is a nearly universal component of the gene regulatory toolkit.
Impact of Confounding Factors	Uncorrected batch effects can inflate false discovery rates (FDR) to ~20.8%, versus an expected 5%.	Technical confounders are pervasive and must be corrected for reliable splicing analysis.

A Complete Workflow for Splicing and Functional Analysis from Bulk RNA-seq

This protocol details the steps from raw sequencing data to the biological interpretation of alternative splicing events, integrating quality control, differential analysis, and functional enrichment.

Protocol: End-to-End Splicing-Aware RNA-seq Analysis

Experimental & Computational Workflow

Step 1: RNA Isolation, Library Preparation, and Sequencing Begin with high-quality total RNA extraction from cells or tissues using a silica-membrane-based kit (e.g., RNeasy Plus) [62]. Assess RNA integrity (RIN > 8 recommended). For mRNA sequencing, use a poly-A selection-based library prep kit (e.g., KAPA mRNA Hyperprep) to generate stranded, paired-end libraries [62]. Sequence on an Illumina platform (HiSeq 4000, NovaSeq) to a minimum depth of 30-40 million paired-end reads per sample for reliable splicing junction detection [15].

Step 2: Read Alignment and Splicing-Aware Quantification Perform initial quality assessment of FASTQ files with FastQC. Align reads to a reference genome (e.g., GRCh38/hg38) using a splice-aware aligner such as STAR or HISAT2 [63]. For differential expression, summarize gene-level counts using featureCounts [62]. For splicing analysis, process the aligned BAM files with a specialized tool:

• Using rMATS: Quantify predefined alternative splicing event types (skipped exon, alternative 5'/3' splice site, mutually exclusive exons, retained intron) by calculating Percent Spliced-In (Ψ, PSI) values for each event [62].
• Using MAJIQ v2: Build a splice graph that includes unannotated junctions and define Local Splicing Variations (LSVs). This approach is particularly suited for large, heterogeneous datasets (hundreds to thousands of samples) and captures complex variation patterns [15].

Step 3: Differential Analysis and Confounder Correction Perform differential expression analysis on gene count matrices using a negative binomial model (DESeq2, edgeR) [62]. For differential splicing, test for significant changes in PSI (ΔΨ) between conditions using the statistical model of the quantification tool (e.g., rMATS uses a likelihood ratio test) [62]. Critically, before differential testing, use the sva R package or similar to identify and correct for batch effects and other unwanted technical variation, which are pervasive and can severely inflate false positive rates if unaddressed [61].

Step 4: Integrated Functional Enrichment with pairedGSEA Translate gene and splicing event lists into biological insight using functional enrichment analysis.

• For differentially expressed genes, perform standard Gene Set Enrichment Analysis (GSEA) using the fgsea package against databases like GO, KEGG, or Reactome.
• For differentially spliced genes, do not rely solely on gene-level enrichment. Use a splicing-aware tool like pairedGSEA [61]. This tool performs GSEA using a statistic derived from the splicing analysis output (e.g., from DEXSeq), allowing direct interrogation of which pathways are enriched for splicing changes.
• Compare and contrast the results from the expression and splicing enrichment analyses. Pathways significantly enriched in both analyses indicate co-regulation, while those unique to one analysis highlight distinct regulatory mechanisms [61].

Step 5: Validation Validate key alternative splicing events using reverse transcription PCR (RT-PCR) or quantitative RT-PCR with primers flanking the alternative exon or junction. Visualize products on a high-resolution gel or via capillary electrophoresis to confirm the predicted isoform shifts [62].

Methodologies for Functional Enrichment of Splicing Data

Moving from lists of significant splicing events to biological understanding requires specialized enrichment methodologies that account for the quantitative nature of splicing changes.

The pairedGSEA Methodology

The pairedGSEA R package is designed for the joint analysis of differential expression and splicing [61].

• Input: It takes results from differential expression (e.g., DESeq2 output) and differential splicing (e.g., DEXSeq output) analyses.
• Ranking: Genes are ranked based on the strength of evidence for differential splicing, independent of their expression change.
• Enrichment Test: It uses a pre-ranked GSEA algorithm (from the fgsea package) to test for the coordinated shift of this splicing-based gene rank within predefined gene sets (pathways) [61].
• Output: The core output is an enrichment score (ES) and normalized enrichment score (NES) for each gene set, indicating whether the genes in that pathway are disproportionately affected by alternative splicing.

Advanced Analysis of Heterogeneous Datasets with MAJIQ HET

For large-scale studies with inherent biological or technical heterogeneity (e.g., GTEx, tumor cohorts), the standard assumption of a uniform PSI per condition breaks down. The MAJIQ HET algorithm within MAJIQ v2 addresses this [15]:

• It first estimates PSI (Ψ) and its confidence interval for each LSV in each sample individually using a Bayesian model.
• It then applies robust, non-parametric rank-based statistical tests (e.g., Mann-Whitney U, InfoScore) to assess significant ΔΨ between sample groups.
• This approach increases power and reproducibility in heterogeneous data by not assuming a single population mean PSI for a condition [15].

Table 3: Comparison of Core Functional Enrichment Strategies for Splicing Data

Method	Core Principle	Ideal Use Case	Key Advantage	Tool/Implementation
Splicing-aware GSEA	Applies pre-ranked GSEA to gene lists ranked by splicing change significance.	Standard bulk RNA-seq experiments with controlled replicates.	Directly links splicing changes to pathways; allows comparison with expression GSEA.	pairedGSEA R package [61]
Over-representation Analysis (ORA)	Tests if genes with significant splicing changes are over-represented in a gene set.	Preliminary, high-level screening of splicing impacts.	Simple, intuitive, and widely available in all enrichment suites.	clusterProfiler, WebGestalt
Module-Based Analysis	Groups complex LSVs into co-regulated "modules" before functional annotation.	Large, discovery-focused studies with many unannotated events.	Reduces complexity, reveals higher-order organization of splicing regulation.	VOILA Modulizer in MAJIQ v2 [15]

Table 4: Research Reagent Solutions for Splicing-Focused RNA-seq

Category	Item	Function & Rationale	Example/Supplier
Wet-Lab Reagents	RNeasy Plus Micro/Mini Kit	Silica-membrane-based total RNA isolation. Removes genomic DNA via a gDNA eliminator column, critical for accurate RNA-seq.	Qiagen [62]
	KAPA mRNA Hyperprep Kit	Library preparation from poly-A-selected RNA. Optimized for Illumina platforms, provides strand specificity.	Roche [62]
	High-Sensitivity DNA/RNA Assay Kits	Accurate quantification and quality assessment of input RNA and final libraries via fluorometry. Essential for sequencing success.	Agilent Bioanalyzer/TapeStation, Qubit (Thermo Fisher)
Core Bioinformatics Tools	STAR Aligner	Splice-aware aligner for RNA-seq reads. Fast, accurate, and widely used for mapping reads across exon junctions.	GitHub: alexdobin/STAR [63]
	rMATS	Detects and quantifies differential alternative splicing from RNA-seq BAM files. User-friendly for standard event types.	RNA-seq MATLAB Toolsuite [62]
	MAJIQ v2 / VOILA v2	Quantifies complex and de novo splicing variations (LSVs). Handles large, heterogeneous datasets and provides advanced visualization.	majiq.biociphers.org [15]
R/Packages for Analysis	pairedGSEA	Performs comparative gene set enrichment analysis for paired differential expression and splicing data.	R package [61]
	DEXSeq	Identifies differential exon usage from RNA-seq count data. A foundational tool for splicing analysis in R.	Bioconductor package [61]
	sva (Surrogate Variable Analysis)	Identifies and estimates surrogate variables for batch effects and other unwanted variation. Critical for confounding factor correction.	Bioconductor package [61]
Reference Databases	GENCODE / Ensembl	High-quality, comprehensive genome annotation. Provides the gene models and transcript structures essential for splicing analysis.	gencodegenes.org, ensembl.org
	MSigDB (Molecular Signatures Database)	Curated collections of gene sets for GSEA (Hallmarks, C2 pathways, C3 regulatory targets, etc.).	GSEA Broad Institute

Interpreting Results: From Pathways to Biological Meaning and Therapeutic Insight

The final, critical step is synthesizing enrichment results into a coherent biological narrative with potential translational relevance.

Conceptual Model for Splicing-Mediated Regulation

Alternative splicing can rewire pathways through several mechanistic consequences that enrichment analysis helps uncover:

• Gain or Loss of Functional Domains: Exon skipping can remove or add protein interaction domains, phosphorylation sites, or localization signals.
• Introduction of Dominant-Negative Isoforms: Switching to a truncated isoform can inhibit the function of the full-length protein from the same allele [61].
• Altered Protein Stability or Localization: Isoforms with different 3' UTRs may have distinct miRNA binding sites, affecting mRNA stability and translation.

Key Interpretation Workflow:

• Identify High-Confidence Pathways: Focus on pathways with strong enrichment scores (FDR < 0.05) from the splicing-aware analysis (e.g., pairedGSEA).
• Characterize Splicing Changes within Pathway: Drill down into the specific genes driving the enrichment. Examine the type of splicing event (e.g., exon skipping, intron retention) and the predicted molecular consequence for each protein (domain loss, PTC introduction).
• Integrate with Expression Data: Determine if the pathway is also transcriptionally regulated. Exclusive splicing regulation suggests precise, post-transcriptional tuning. Co-regulation suggests a multi-layered response.
• Formulate a Testable Hypothesis: Synthesize findings into a model. For example: "Treatment X induces a splicing shift in the apoptosis pathway, favoring pro-apoptotic isoforms of BCL2L1 and CASP2, which coincides with observed cell death, independent of transcriptional changes in these genes."
• Prioritize for Therapeutic Development: In drug discovery, pathways uniquely regulated by splicing in disease tissues represent potential novel targets for small molecules modulating splice factors or antisense oligonucleotides (ASOs) designed to correct specific aberrant splicing events.

Disease Context and Transcriptional Complexity

The analysis of alternative splicing (AS) via bulk RNA sequencing (RNA-seq) provides a critical window into the molecular pathology of complex human diseases. In cancer and neurological disorders, widespread dysregulation of splicing alters proteome diversity, drives oncogenic transformation, and contributes to neuronal dysfunction. This post-transcriptional regulatory layer, affecting over 90% of multi-exon human genes, serves as a rich source of diagnostic biomarkers and therapeutic targets [8]. Framed within a broader thesis on bulk RNA-seq research, this document details application notes and protocols for investigating splicing dysregulation in these disease contexts, emphasizing reproducible methodologies for researchers and drug development professionals.

1.1 Splicing Dysregulation in Claudin-Low Breast Cancer Triple-negative breast cancers (TNBC) are aggressive malignancies with poor prognosis. The Claudin-low subtype, constituting 25–39% of TNBC cases, is defined by low cell adhesion protein expression and mesenchymal features [64]. A pivotal study investigating the transmembrane co-receptor Neuropilin-1 (NRP1) in Claudin-low breast cancer provides a model for splicing-centric oncology research. Knockdown of NRP1 in three cell lines (HS578T, MDA-MB-231, SUM159PT) followed by bulk RNA-seq (52 samples) generated a dataset for transcriptomic and isoform-level analysis (GEO: GSE266566) [64]. NRP1 is implicated in regulating tumor progression through stem cell characteristics and RAS/MAPK pathway activation [64]. Analyzing this dataset for AS events—such as exon skipping in genes related to epithelial-mesenchymal transition (EMT) or therapy resistance—can reveal how splicing factors and signaling pathways converge to drive an aggressive cancer phenotype.

1.2 Splicing Alterations in Alzheimer’s Disease In Alzheimer’s disease (AD), transcriptomic changes precede and accompany pathological hallmarks like amyloid plaques. A systematic meta-analysis of bulk RNA-seq studies on human AD brain tissue identified 571 differentially expressed genes (DEGs) in the temporal lobe and 189 in the frontal lobe compared to non-demented controls [65]. Key pathways, including "neuroactive ligand-receptor interaction," were highlighted [65]. Independent analysis of 221 samples (132 AD, 89 control) identified 12 robust DEGs, with Transthyretin (TTR) showing significant downregulation and a strong link to amyloid fiber formation pathways [66]. Beyond differential expression, alternative splicing events likely affect genes central to neuronal function, synaptic integrity, and immune response in the brain. For instance, altered splicing of pre-synaptic genes or tau (MAPT) could directly contribute to pathophysiology. Investigating splicing patterns in brain region-specific compendia, such as those harmonized from 50 sources encompassing 16,446 datasets, is essential for distinguishing driver splicing events from bystander effects [67].

Key Experimental and Computational Protocols

2.1 Protocol 1: Bulk RNA-seq for Splicing Analysis from Cell Line Models This protocol outlines the generation of a splicing-focused RNA-seq dataset, based on the Claudin-low breast cancer study [64].

• Step 1: Experimental Perturbation and RNA Extraction.

  • Culture relevant cell lines (e.g., HS578T, MDA-MB-231) under standard conditions.
  • Perform knockdown of the gene of interest (e.g., NRP1) using two distinct siRNA or shRNA sequences alongside a non-targeting (NT) control. Include biological replicates (minimum n=3).
  • After 48-72 hours, extract total RNA using a column-based or magnetic bead purification kit. Assess RNA integrity (RIN > 8.0) using an Agilent Bioanalyzer or TapeStation [67].

• Step 2: Library Preparation and Sequencing.

  • Use 50–100 ng of high-quality total RNA as input [67].
  • Perform ribosomal RNA depletion using enzymatic methods (e.g., Illumina Ribo-Zero Plus) to retain non-coding RNAs and unprocessed transcripts, which can inform on intron retention events [67].
  • Construct stranded RNA-seq libraries using a ligation-based kit (e.g., Illumina Stranded Total RNA Prep). This preserves strand information, crucial for accurate splice junction assignment.
  • Perform paired-end sequencing (100 bp x 2) on an Illumina platform to a minimum depth of 30-40 million mapped reads per sample for reliable splicing quantification [64] [67].

• Step 3: Core Bioinformatics Processing.

  • Quality Control: Use FastQC to assess raw read quality. Trim adapters and low-quality bases with Trimmomatic or Cutadapt.
  • Alignment: Align reads to the human reference genome (GRCh38) using a splice-aware aligner such as STAR [67].
  • Quantification: Generate two levels of quantification:
    • Gene-level counts: Using featureCounts or HTSeq for differential expression analysis [66].
    • Transcript/Isoform-level counts: Using a tool like RSEM (as used in the Treehouse compendia) or Salmon for isoform abundance estimation [67].
  • Splice Junction Analysis: Extract novel and annotated splice junctions directly from the STAR alignment files for downstream event detection.

2.2 Protocol 2: Computational Identification of Differential Splicing Events This protocol details the analysis of processed RNA-seq data to identify statistically significant alternative splicing events.

• Step 1: Event Detection and Quantification.

  • Use specialized software to identify and quantify specific AS event types (exon skipping, intron retention, alternative 5'/3' splice sites).
  • Recommended Tools: rMATS, SUPPA2, or LeafCutter. These tools calculate Percent Spliced In (PSI or Ψ) values for each event in every sample.
  • Input: Aligned BAM files or junction read counts, along with a reference transcriptome annotation (GTF file).

• Step 2: Differential Splicing Analysis.

  • Perform statistical testing to compare PSI values between experimental groups (e.g., NRP1-KD vs. NT control; AD vs. control brain tissue).
  • Tools like rMATS employ likelihood-ratio tests to identify differential splicing events with a false discovery rate (FDR) correction. A common threshold is FDR < 0.05 and |ΔPSI| > 0.1.
  • For Alzheimer's meta-analysis, harmonize results across multiple studies using random-effects models, as demonstrated in [65].

• Step 3: Integration with Differential Expression and Pathway Analysis.

  • Integrate lists of differentially spliced genes (DSGs) with differentially expressed genes (DEGs) from the same dataset.
  • Perform functional enrichment analysis (Gene Ontology, KEGG) on the combined or unique gene sets using tools like Enrichr [66]. This identifies biological pathways most impacted by transcriptional and post-transcriptional dysregulation (e.g., "RAS/MAPK signaling" in cancer, "amyloid fiber formation" in AD).
  • Visualize networks of interacting proteins encoded by DSGs/DEGs using STRING and Cytoscape to identify hub genes [66].

2.3 Protocol 3: Analysis of Alternative Polyadenylation (APA) APA is a coupled regulatory mechanism that influences mRNA stability, localization, and protein diversity. Its analysis complements AS studies [2].

• Step 1: Data Processing for APA.

  • The input can be standard bulk RNA-seq data. The deep-learning tool PolyAMiner-Bulk is specifically designed for this purpose [68].
  • Follow the tool's pre-processing requirements, which typically involve using a specific aligner (e.g., STAR) and generating certain file formats (BAM, genome coverage files).

• Step 2: Running PolyAMiner-Bulk.

  • Install the tool in a high-performance computing environment with GPU support recommended [68].
  • Execute the pipeline, which uses a deep neural network to decode APA dynamics, identifying and quantifying proximal and distal poly(A) site usage.
  • The output includes metrics like the Relative Expression Difference (RED) for APA events between conditions.

• Step 3: Downstream Interpretation.

  • Generate result visualizations such as APA landscape plots, PCA plots, and volcano plots to assess global APA changes [68].
  • Integrate significant APA events with DSGs and DEGs to build a comprehensive view of post-transcriptional deregulation in the disease model.

Data Presentation and Synthesis

3.1 Quantitative Summary of Splicing Dysregulation Table 1: Key Quantitative Findings from Featured Studies

Disease Context	Dataset/Source	Key Splicing/Expression Metric	Result	Implication
Claudin-Low Breast Cancer	Bulk RNA-seq of NRP1 knockdown (GSE266566) [64]	Differentially expressed genes (DEGs) & isoforms	Comprehensive transcriptomic landscape for 3 cell lines post-NRP1 KD.	Identifies NRP1-regulated pathways (e.g., RAS/MAPK) as potential therapeutic targets.
Alzheimer’s Disease (Meta-Analysis)	Systematic review of 24 studies [65]	Number of significant DEGs (FDR < 0.05)	571 DEGs in temporal lobe; 189 DEGs in frontal lobe.	Highlights brain region-specific vulnerability and potential biomarkers (e.g., CLU, GFAP).
Alzheimer’s Disease (Biomarker Study)	Analysis of 221 patient samples [66]	Number of prioritized DEGs & top candidate	12 robust DEGs identified; TTR most downregulated.	Suggests TTR as a diagnostic biomarker and links it to amyloid pathology.
Pan-Cancer Pediatric Resource	Treehouse Compendia [67]	Number of harmonized RNA-seq datasets	16,446 datasets, 5,687 from pediatric/adolescent/young adult individuals.	Enables powerful comparative splicing analysis for rare cancers.

3.2 The Scientist's Toolkit: Essential Reagents and Resources Table 2: Key Research Reagent Solutions for Splicing Analysis

Item Name	Function/Application	Example/Provider
Ribo-Zero Plus rRNA Depletion Kit	Removes ribosomal RNA from total RNA samples, preserving pre-mRNA and non-coding RNA for full transcriptome analysis.	Illumina [67]
Stranded Total RNA Prep, Ligation with Ribo-Zero Plus Kit	Library preparation kit for strand-specific, rRNA-depleted RNA-seq libraries. Essential for accurate splice junction orientation.	Illumina [67]
siRNA/shRNA for Gene Knockdown	Induces targeted gene silencing (e.g., NRP1) to study its effect on the global transcriptome and splicing landscape.	Various (e.g., Dharmacon) [64]
PolyAMiner-Bulk Software	Deep-learning-based algorithm to identify and quantify Alternative Polyadenylation (APA) events from bulk RNA-seq data.	[68]
Toil RNA-seq Pipeline	Dockerized, portable workflow for consistent processing of RNA-seq data, including alignment (STAR) and quantification (RSEM).	UCSC Treehouse [67]
Longcell Computational Pipeline	Specialized tool for isoform quantification and splicing analysis from single-cell/spatial Nanopore long-read data.	[8]
Harmonized Expression Compendia	Pre-processed, uniformly analyzed bulk RNA-seq datasets from thousands of samples for comparative meta-analysis.	UCSC Treehouse (5 compendia) [67]

Signaling Pathways and Therapeutic Implications

4.1 NRP1 Signaling in Cancer Progression Neuropilin-1 (NRP1) is a pleiotropic co-receptor that enhances signaling of pathways like VEGF, TGF-β, and EGFR [64]. In Claudin-low breast cancer, NRP1 overexpression promotes tumor progression through acquisition of stem cell characteristics and RAS/MAPK pathway activation [64]. Splicing analysis can reveal how NRP1 knockdown alters the isoform balance of key components within these pathways, potentially modulating their activity. For example, alternative splicing of genes in the MAPK cascade can create isoforms with altered kinase activity or subcellular localization, adding a layer of regulation to this oncogenic pathway.

4.2 From Splicing Insights to Therapeutic Strategies Dysregulated splicing generates unique disease-associated protein isoforms, which represent attractive therapeutic targets.

• Small Molecules & Antisense Oligonucleotides (ASOs): Identified splicing events can be targeted directly. For example, an exon-skipping event creating a pro-oncogenic isoform could be corrected using ASOs. In neurological disorders, ASOs are already in clinical trials for conditions like spinal muscular atrophy and are a promising avenue for AD.
• RNA Cancer Vaccines: Neoantigens derived from cancer-specific splice variants are ideal targets for personalized immunotherapy. Recent breakthroughs in RNA vaccine technology, including platforms using circular RNA and novel lipid nanoparticles, show efficacy in clinical trials for melanoma and pancreatic cancer [69]. The manufacturing timeline for personalized vaccines has been reduced to under four weeks [69]. Splicing analysis of patient tumors can identify the most immunogenic, tumor-specific isoforms for vaccine design.
• Drug Repurposing: Computational analysis linking splicing signatures to drug response can identify repurposing opportunities. In the AD study, molecular docking suggested the FDA-approved drug Levothyroxine could interact with the downregulated TTR protein [66].

Visualizing Workflows and Pathways

Bulk RNA-Seq Splicing Analysis Workflow

NRP1 Signaling and Splicing Crosstalk in Cancer

Optimizing Your Splicing Analysis: Overcoming Challenges for Robust Results

Empirical Guidelines for Bulk RNA-Seq Splicing Analysis

Designing a robust bulk RNA-seq experiment for alternative splicing analysis requires balancing statistical power, sensitivity to detect isoform-level changes, and practical constraints. Empirical data and simulation studies provide concrete guidance on the interplay between sample size, sequencing depth, and analytical choice.

Quantitative Trade-offs: Replicates, Depth, and Detection Power

A foundational simulation study based on six diverse public RNA-seq datasets provides critical benchmarks for power [70]. The analysis demonstrates that increasing the number of biological replicates confers a greater boost to statistical power than increasing sequencing depth per sample, especially once a depth of approximately 20 million reads is achieved [70]. For a fixed budget, a design favoring more samples at a moderate depth often yields locally optimal power [70].

The performance of differential expression (DE) analysis software varies. In comparative evaluations, DESeq2 and edgeR tend to provide the best overall performance for gene-level analysis [70]. For experiments with very small sample sizes, the sSeq package may offer better sensitivity, while EBSeq is noted for its robustness in isoform-level expression analysis [70]. It is critical to note that detection power is lower for long intergenic noncoding RNAs (lincRNAs) compared to protein-coding mRNAs due to their generally lower expression levels [70].

Table 1: Comparative Performance of RNA-Seq Differential Expression Analysis Packages [70]

Analysis Package	Recommended Use Case	Key Performance Note
DESeq2	General gene-level DE; Paired designs	High overall power and performance.
edgeR	General gene-level DE	High overall power and performance.
EBSeq	Isoform-level expression analysis	Robust for isoforms; not adapted for paired-sample designs.
sSeq	Experiments with very small sample sizes	Better sensitivity for low replication.

For the specific challenge of detecting alternative splicing events, the required sample size escalates. A probabilistic model for observing splice junctions (gap-sites) indicates that events are heterogeneously distributed, with many rare, sporadically used splice sites [71]. Consequently, merging data from multiple samples is essential to capture this diversity. In a study of 54 human fibroblast transcriptomes, the number of detected high-confidence splice sites (gap-quality level 3) grew slowly as sample batches increased from 2 to 12, while low-confidence events increased nearly linearly [71]. This underscores that detecting a comprehensive repertoire of splicing events, especially rarer ones, demands larger sample cohorts.

Table 2: Sample Size and Splicing Event Detection (Based on Fibroblast Data) [71]

Number of Samples in Batch	Total Gap-Sites Detected (approx.)	High-Confidence (gql3) Gap-Sites (approx.)
2	706,000	92,000
4	1,076,000	105,000
8	1,708,000	124,000
12	2,270,000	137,000

Core Statistical Foundations: Error Rates and Effect Size

Sample size calculation is intrinsically linked to managing statistical error rates and the anticipated biological effect [72].

• Type I Error (α): The probability of falsely rejecting a true null hypothesis (false positive). A standard threshold is α=0.05 [72].
• Type II Error (β): The probability of failing to reject a false null hypothesis (false negative) [72].
• Statistical Power (1-β): The probability of correctly rejecting a false null hypothesis. A power of 0.8 (80%) is conventionally considered adequate [72] [73].
• Effect Size (ES): The magnitude of the biological difference (e.g., fold-change in expression or percent spliced in). A smaller effect size requires a larger sample size to detect it with the same power [72] [73].

The relationship between these factors is delicate: reducing the Type I error rate (α) or increasing power (1-β) necessitates a larger sample size, as does aiming to detect a smaller effect size [73].

Protocols for Sample Size Calculation and Experimental Design

Protocol: Sample Size Calculation for Differential Splicing

This protocol outlines a method for calculating sample size per group for a bulk RNA-seq experiment aiming to detect differential splicing events, based on a negative binomial model and controlling for the False Discovery Rate (FDR) [74].

Step 1: Define Statistical Parameters.

• Set the significance level (α), commonly 0.01 or 0.05 [72].
• Define the desired statistical power (1-β), typically 0.8 or 0.9 [72] [74].
• Establish the target False Discovery Rate (FDR) for multiple testing correction (e.g., 0.05) [74].

Step 2: Estimate Biological Parameters from Pilot Data or Literature.

• Mean expression (μ₀): Estimate the average normalized read count (e.g., in Counts Per Million) for the splicing event or gene of interest in the control group [74].
• Dispersion (φ): Estimate the biological variance of the event. This can be derived from a pilot study or public datasets with similar biological material [70] [74].
• Fold Change (ρ): Define the minimum alternative splicing ratio or gene expression fold change you wish to detect (e.g., 1.5, 2.0) [74].
• Normalization factor ratio (w): Account for differences in library size between groups [74].

Step 3: Employ a Computational Calculation Tool. Due to the complexity of the negative binomial model, analytical solutions are impractical. Use one of the following:

• R-based functions: Implement the power calculation based on the exact test for negative binomial data, as described in [74]. This involves iterative numerical methods to solve for sample size n in the power function 1-β = ξ(n, ρ, μ₀, φ, w, α).
• Online calculators: Use tools like the RNA-Seq Power Calculator (based on [70]) or other statistical power calculators that incorporate dispersion estimates [70] [73].
• Simulation: For complex designs, perform a Monte Carlo simulation: simulate read count data under the alternative hypothesis using parameters from Step 2, analyze it with your chosen pipeline (e.g., edgeR, DEXSeq), and calculate the proportion of simulations where the event is correctly called significant. Adjust sample size until the proportion meets your target power [70].

Step 4: Adjust for Multiple Testing and Splicing Complexity.

• The sample size calculated in Step 3 is for a single event. To control the FDR across thousands of tested splicing events, the required sample size will be larger [74].
• For splicing analysis, consider that many events are rare. To model the probability of observing a splicing event in at least k samples, a binomial model can be applied: P(observation in ≥k samples) = 1 - Σ_{i=0}^{k-1} [(n choose i) * p^i * (1-p)^{n-i}], where p is the per-sample observation probability of that event [71]. Use this to determine n needed to observe rare events with high probability.

Protocol: Optimal Experimental Design for Splicing Studies

A. Replication Strategy

• Biological Replicates: These are non-negotiable. They capture random biological variation and are the primary unit for statistical analysis [75] [31]. A minimum of 3 biological replicates per condition is required, but 4-8 is strongly recommended for reliable inference, especially for heterogeneous samples or to detect small effect sizes [31] [76].
• Technical Replicates: With modern RNA-seq protocols, technical variation is low. Technical replicates are generally unnecessary for quantifying biological differences and are a poor use of resources [75].

B. Sequencing Depth Guidance

• General Gene-level DE: 15-30 million aligned reads per sample is often sufficient when paired with an adequate number of replicates [75] [76].
• Isoform/Splicing Analysis: Deeper sequencing is required. A minimum of 30 million paired-end reads per sample is recommended for known isoforms, and >60 million for novel isoform discovery [75]. Longer read lengths (≥100 bp paired-end) improve the mapping of reads across exon junctions [75] [76].

C. Avoiding Confounding and Batch Effects

• Confounding: Ensure biological variables (e.g., sex, age, batch) are not perfectly correlated with your experimental condition. Balance these factors across treatment groups [75].
• Batch Effects: Systematic technical differences between processing batches are a major source of false findings [75].
  • Design: Distribute samples from all experimental conditions across every processing batch (RNA extraction, library prep, sequencing run) [75].
  • Documentation: Meticulously record all batch information in sample metadata [75].
  • Analysis: Use statistical models (e.g., in DESeq2) to include "batch" as a covariate during differential testing to regress out its effect [75].

Workflow Visualization for Study Design and Analysis

Diagram 1: Sample Size Determination & Experimental Design Workflow (97 characters)

Diagram 2: Bulk RNA-Seq Splicing Analysis Pipeline (86 characters)

The Scientist's Toolkit: Essential Reagents and Materials

Table 3: Research Reagent Solutions for RNA-Seq Splicing Studies

Item / Solution	Function / Purpose	Key Considerations
High-Quality Total RNA	Starting material for library preparation.	RNA Integrity Number (RIN) >7 is critical for mRNA-Seq; RIN >5 may be acceptable for total RNA-Seq protocols [76].
Poly-dT Enrichment or rRNA Depletion Kits	Selection of RNA species.	Poly-dT enrichment is standard for poly-A mRNA. rRNA depletion is required for total RNA analysis, including non-coding RNAs, and is more tolerant of moderate RNA degradation [76].
Stranded Library Prep Kit	Creates sequencing libraries that preserve the original strand orientation of transcripts.	Essential for accurate annotation of transcripts and identification of antisense transcription [76].
External RNA Spike-in Controls (e.g., SIRVs)	Added to samples in known concentrations to monitor technical performance.	Enables assessment of dynamic range, sensitivity, and quantification accuracy across the entire workflow [31].
DNase I	Removal of contaminating genomic DNA.	A critical step during RNA purification or as part of the library prep protocol to prevent false-positive signals [76].
Dual-Indexed Adapters	Unique molecular identifiers for each sample.	Allows for multiplexing of many samples in a single sequencing lane and correct demultiplexing, essential for cost-effective large-scale studies [76].
Tumor Tissue Microarray (TMA)	Validation platform for splicing biomarkers.	Enables high-throughput validation of candidate splicing events or protein isoforms via RNA in situ hybridization (RNA-ISH) or immunohistochemistry across hundreds of tissue specimens simultaneously [77].
Automatic Tissue Microarrayer	Construction of TMAs for validation studies.	Precisely extracts tissue cores from donor blocks and assembles them into recipient blocks, standardizing validation workflows [77] [78].

Handling Technical and Biological Heterogeneity in Large Datasets

Within the broader thesis context of alternative splicing analysis using bulk RNA sequencing, addressing both technical and biological heterogeneity is a fundamental challenge. Large-scale transcriptomic studies, such as those involving thousands of samples from consortia like GTEx, inherently contain multiple sources of variation [15]. Technical heterogeneity arises from differences in sample preparation, library protocols, sequencing platforms, and batch effects [19]. Biological heterogeneity stems from the diversity of cell types within a tissue, individual genetic variation, and dynamic biological states [79]. In splicing analysis, this heterogeneity manifests as increased variability in the Percent Spliced In (PSI, Ψ) values for local splicing variations (LSVs), complicating the detection of true differential splicing events [15]. Traditional bulk RNA-seq provides an average gene expression profile across thousands to millions of cells, obscuring cell-type-specific splicing patterns [80] [81]. Emerging computational and experimental strategies are now focused on deconvolving these complex signals, improving the accuracy and biological relevance of splicing analyses in large, heterogeneous datasets.

Technical noise is introduced at every stage of an RNA-seq experiment. Key sources include:

• Library Preparation: Variations in reverse transcription efficiency, amplification bias (PCR or IVT), and adapter ligation [82] [81].
• Sequencing: Differences in sequencing depth, platform-specific biases, and run-to-run variability [19].
• Sample Processing: Batch effects from different researchers, isolation days, or reagent lots. Inadequate quality control of input RNA (RIN, concentration) also contributes significantly [19].

Biological heterogeneity presents a more complex layer:

• Cellular Composition: Bulk tissue RNA-seq averages expression across all constituent cell types. Changes in cell type proportions between conditions can be misidentified as changes in gene expression or splicing within a cell type [80].
• Individual Variation: Genetic background, age, sex, and environmental exposures lead to expression and splicing differences between individuals [15].
• Stochastic Gene Expression: Even within a homogeneous cell population, intrinsic noise leads to cell-to-cell variation in transcript abundance [79].

Impact on Splicing Analysis

For alternative splicing analysis, heterogeneity is particularly problematic. Splicing metrics like PSI are ratios, making them sensitive to the underlying expression of the gene and the specific splice junctions [15]. In large, heterogeneous datasets, the assumption that samples within a condition share a common PSI value is often violated. This can lead to:

• Loss of Statistical Power: Increased within-group variance reduces the ability to detect true differential splicing.
• Increased False Positives: Uncorrected batch effects or confounding biological variables can create spurious signals of splicing change [15].
• Reduced Reproducibility: Findings may not replicate across studies with different technical or population structures.

Computational Strategies and Tools for Managing Heterogeneity

A suite of computational tools has been developed to model, quantify, and correct for heterogeneity in large RNA-seq datasets. The choice of tool depends on the analysis goal (e.g., differential expression vs. splicing) and the study design.

Table 1: Computational Tools for Addressing Heterogeneity in RNA-Seq Analysis

Tool/Method	Primary Purpose	Approach to Heterogeneity	Key Features	Reference
MAJIQ v2 (HET)	Differential Splicing	Non-parametric, rank-based tests (TNOM, Mann-Whitney U) that do not assume a common PSI per group. Quantifies PSI per sample before group comparison.	Handles large, heterogeneous datasets; incremental splicegraph builder; classifies complex LSVs into modules.	[15]
DiSC	Differential Expression (scRNA-seq)	Models individual-level biological variability by extracting multiple distributional characteristics and using a permutation testing framework.	Scalable to large sample sizes; ~100x faster than some state-of-the-art methods; controls FDR.	[83]
gCCA (genoMap CCA)	Cellular Deconvolution	Uses a deep learning framework (convolutional VAE) on 2D gene-gene interaction images to derive sample-specific cell type signatures, reducing noise susceptibility.	Accounts for inter-sample variability; improves deconvolution accuracy by >14% on average.	[80]
CIBERSORTx	Cellular Deconvolution	Uses support vector regression (SVR) with predefined gene signature matrices to infer cell type abundances from bulk data.	A widely used classical method; requires a pre-defined signature matrix.	[80]
Curare/GenExVis	Workflow & Visualization	A modular, reproducible Snakemake-based workflow builder that standardizes preprocessing, alignment, and DGE analysis to minimize technical variability in processing.	Ensures reproducible analysis; customizable modules for various aligners and tools; coupled with an offline visualization GUI.	[84]
Salmon	Expression Quantification	Uses "pseudoalignment" and a robust statistical model to account for uncertainty in read assignment to transcripts, crucial for accurate isoform-level quantification.	Fast, alignment-free mode available; can also use alignment files (BAM); integrates well with workflows.	[18]

Workflow for Heterogeneity-Aware Splicing Analysis

The following diagram outlines a recommended computational workflow for analyzing alternative splicing in large, heterogeneous bulk RNA-seq datasets, incorporating strategies to mitigate technical and biological noise.

Diagram 1: Workflow for Splicing Analysis in Heterogeneous Datasets

Statistical Framework for Heterogeneous Splicing

MAJIQ v2 introduces a new statistical paradigm (MAJIQ HET) to handle heterogeneous groups. The following diagram contrasts the traditional homogeneous group model with the new heterogeneous sample model.

Diagram 2: Statistical Models for Splicing in Homogeneous vs. Heterogeneous Groups

Detailed Experimental Protocols

Protocol: Bulk RNA-Seq with Sample Multiplexing for Batch Mitigation

This protocol minimizes technical batch effects by using sample multiplexing (barcoding) early in the workflow, allowing multiple samples to be processed as a single library [81].

Objective: Generate high-quality, batch-minimized sequencing libraries from multiple bulk RNA samples for downstream splicing analysis. Materials: See "The Scientist's Toolkit" (Section 6). Procedure:

• RNA Extraction & QC: Extract total RNA using a column-based or TRIzol method. Assess RNA integrity (RIN > 8.0 recommended) and concentration using a TapeStation/Bioanalyzer and Qubit [81] [19].
• Sample Normalization: Dilute all samples to the same RNA concentration (e.g., 10 ng/μL) to minimize read-depth variability [81].
• Reverse Transcription with Unique Barcodes: For each sample, perform first-strand cDNA synthesis using a primer containing a unique sample barcode. Pool samples after this step for all subsequent reactions [82] [81].
  • Critical Step: Use Unique Molecular Identifiers (UMIs) within the RT primer to correct for amplification bias in downstream quantification [82] [79].
• cDNA Amplification & Library Prep: Amplify the pooled, barcoded cDNA using a high-fidelity polymerase. Proceed with standard library preparation: fragment cDNA, perform end-repair/A-tailing, and ligate sequencing adapters.
• Quality Control & Sequencing: Validate the final library size distribution and concentration (e.g., via TapeStation). Sequence on an Illumina platform using paired-end reads (recommended: 2x100bp for junction coverage) [18].

Protocol: Computational Pipeline for Splicing Analysis with MAJIQ v2

This protocol details the analysis of aligned RNA-seq data to detect differential splicing in heterogeneous conditions [15].

Objective: Identify and visualize differential Local Splicing Variations (LSVs) from large, heterogeneous bulk RNA-seq datasets. Input: Aligned RNA-seq reads in BAM format (generated by a splice-aware aligner like STAR). Software: MAJIQ v2 package (including VOILA v2 for visualization). Procedure:

• Build the Splicing Graph:
  • The configuration file (conf.txt) lists BAM files and sample groups. The builder constructs a splice graph for each gene, incorporating annotated and de novo (unannotated) junctions present in the data. Low-coverage junctions are filtered.
• Quantify Splicing (PSI):
  • The quantifier calculates the posterior distribution of PSI (Ψ) for every junction in each LSV for every sample. Use the --het flag to enable sample-level PSI estimation for heterogeneous group analysis.
• Differential Splicing Analysis:
  • For heterogeneous groups, ensure the configuration specifies the use of the HET (TNOM or Mann-Whitney) test statistic. This tests for significant changes in the distribution of sample-level PSI values between groups.
• Module Classification & Visualization:
  • Use the modulize command to classify complex LSVs into coherent alternative splicing modules (e.g., cassette exons, mutually exclusive exons). Use voila view to interactively visualize splice graphs, PSI distributions across hundreds of samples, and significant deltapsi events.

Protocol: Cellular Deconvolution Using genoMap-Based CCA (gCCA)

This protocol estimates cell type proportions from bulk RNA-seq data using a deep learning approach that accounts for inter-sample heterogeneity, providing a biological covariate for splicing analysis [80].

Objective: Estimate sample-specific cell type abundances from bulk RNA-seq data to control for cellular composition effects. Input: A gene expression count matrix from bulk RNA-seq data and a reference scRNA-seq dataset for the tissue of interest. Software: gCCA Python implementation. Procedure:

• Prepare Reference Data: Process the scRNA-seq reference to obtain a normalized gene expression matrix. Identify major cell types via clustering.
• Construct genoMaps: Transform the scRNA-seq expression profiles for each cell into 2D genoMap images, where gene-gene correlations are encoded as spatial proximity.
• Train Convolutional VAE & GMM: Train a convolutional Variational Autoencoder (VAE) on all reference genoMaps. Apply a Gaussian Mixture Model (GMM) in the latent feature space to define sample-specific, cell-type-signature genoMaps.
• Deconvolve Bulk Samples: For each bulk sample, construct its genoMap. In the image domain, decompose this bulk genoMap into a linear combination of the learned cell-type-signature genoMaps to estimate cellular proportions.
• Integration with Splicing Analysis: Use the estimated cell type proportions as covariates in the MAJIQ HET differential splicing model (e.g., included in the sample grouping/design matrix) to control for confounding due to cellular heterogeneity.

Data Management and Quality Control

Effective management of heterogeneity begins with rigorous experimental design and quality control.

Table 2: Key Quality Control Checkpoints and Metrics

Stage	Checkpoint	Metric/Tool	Action Threshold
Wet Lab	RNA Integrity	RIN (TapeStation)	RIN > 8.0 for standard protocols [19].
	Library Complexity	Fragment Analyzer	Sharp peak at expected insert size.
Sequencing	Raw Read Quality	FastQC, MultiQC	Per base Q-score > 30; adapter contamination minimal [84].
Alignment	Mapping Efficiency	STAR Logs	>70% uniquely mapped reads for typical eukaryotic samples [18].
	Splicing Awareness	Junction Saturation	Ensure sufficient read depth for junction detection.
Quantification	Sample Correlation	PCA Plot	Identify outliers; biological groups should cluster [19].
	Batch Effect	PCA/Batch-aware plots	Look for clustering by processing date or lane.
Splicing-Specific	Junction Coverage	MAJIQ Builder Output	Filter junctions with < 10 reads in < 20% of samples per group [15].

Mitigating Batch Effects Proactively: The most effective strategy is experimental:

• Randomization: Process samples from different experimental groups simultaneously.
• Blocking: If processing in batches is unavoidable, ensure each batch contains a similar proportion of samples from all conditions.
• Positive Controls: Use spike-in RNAs or external RNA controls consortium (ERCC) standards to monitor technical performance across batches [19].

The Scientist's Toolkit: Research Reagent Solutions

Table 3: Essential Reagents and Materials for Heterogeneity-Aware RNA-Seq Studies

Item	Function & Relevance to Heterogeneity	Example Product/Kit
UMI-Adapter RT Primers	Contains Unique Molecular Identifiers to label each original mRNA molecule, correcting for amplification bias and improving quantification accuracy [82] [79].	CEL-Seq2 primers; Commercial UMI kits (e.g., from Illumina, Takara).
Sample Multiplexing Barcodes	Unique nucleotide tags for each sample, enabling early pooling to minimize post-pooling technical variability [82] [81].	Chromium i7 Multiplex Kit; IDT for Illumina RNA UD Indexes.
High-Fidelity Polymerase	For cDNA amplification with low error rates, reducing sequence-based noise.	KAPA HiFi HotStart ReadyMix; Q5 High-Fidelity DNA Polymerase.
RNA Integrity Assay Kits	To rigorously assess input RNA quality, a major source of technical variation [19].	Agilent RNA 6000 Nano/Pico Kit (for TapeStation/Bioanalyzer).
RNA Spike-In Controls	Exogenous RNA added at known concentrations to monitor technical variation in library prep and sequencing between samples [19].	ERCC RNA Spike-In Mix (Thermo Fisher).
Single-Cell RNA-Seq Kit	To generate the high-quality reference profile needed for cellular deconvolution protocols [80].	10x Genomics Chromium Single Cell 3' Kit; SMART-Seq HT Kit.
Splice-Aware Aligner Software	Accurately maps reads across exon junctions, foundational for all splicing analysis.	STAR [18]; HISAT2.
Statistical Computing Environment	Platform for implementing differential expression and splicing analysis tools.	R Statistical Environment (with Bioconductor packages); Python.

In high-throughput biology, researchers routinely test thousands of hypotheses simultaneously. In a standard bulk RNA sequencing (RNA-seq) experiment for differential gene expression, this involves a statistical test for each of approximately 20,000 genes [85]. Even when no true biological differences exist, using a conventional significance threshold (α=0.05) would be expected to yield around 1,000 false positive findings by chance alone [86]. This problem is magnified in alternative splicing (AS) analysis, where the unit of testing shifts from genes to individual splicing events or junctions, often increasing the number of hypotheses by an order of magnitude [2].

The core challenge is controlling the inflation of Type I errors (false positives) while maintaining power to detect true biological signals. Traditional corrections like the Bonferroni method control the Family-Wise Error Rate (FWER)—the probability of at least one false positive—but are often deemed excessively conservative for exploratory omics studies, leading to many missed findings (false negatives) [87] [88]. The False Discovery Rate (FDR) has emerged as a more suitable alternative for large-scale studies [87]. Formally, the FDR is defined as the expected proportion of false discoveries among all discoveries made (FDR = E[V/R], where V is the number of false positives and R is the total number of significant findings) [87].

This article details the application of FDR control within the framework of alternative splicing analysis in bulk RNA-seq. AS is a critical post-transcriptional regulator, with over 95% of multi-exon human genes undergoing splicing to produce diverse transcript isoforms [2]. Accurately identifying differential splicing events (DSEs) between conditions is therefore paramount. We outline statistical considerations, provide specific protocols for analysis, and highlight how feature dependency—a hallmark of splicing data—impacts FDR control and interpretation.

Foundational Concepts and Statistical Definitions

Understanding the landscape of error control requires distinguishing between key related metrics, each with a specific interpretation and use case.

Table 1: Key Statistical Metrics for Multiple Testing Correction

Metric	Full Name	Definition	Primary Use Case	Control in Splicing Context
P-value	Probability value	Probability of observing a test statistic as extreme as, or more extreme than, the one calculated, assuming the null hypothesis is true [89].	Single hypothesis test.	Raw significance for an individual splicing event.
PCER	Per-Comparison Error Rate	Expected proportion of false positives out of all hypothesis tests conducted [87].	Uncorrected multiple testing.	Not controlled; leads to rampant false positives.
FWER	Family-Wise Error Rate	Probability of making one or more false discoveries among all hypotheses tested [87].	Confirmatory studies where any false positive is costly.	Controlled by Bonferroni; often too strict for splicing discovery.
FDR	False Discovery Rate	Expected proportion of false discoveries among all discoveries declared significant [90] [87].	Exploratory high-throughput studies (e.g., RNA-seq, AS).	Target for methods like Benjamini-Hochberg; balances discovery and error.
Q-value	-	The minimum FDR at which a given test (e.g., a splicing event) may be called significant [87] [89].	Ranking and thresholding findings based on FDR.	Direct FDR estimate for each event; output of many splicing tools.

The choice between FWER and FDR is a trade-off between stringency and power. The Bonferroni correction controls FWER by dividing the target significance threshold (α) by the number of tests (m): α_adjusted = α / m. For 100,000 splicing events at α=0.05, the adjusted threshold becomes 5 x 10⁻⁷, making detection exceedingly difficult [88]. In contrast, FDR-controlling procedures like Benjamini-Hochberg (BH) are less conservative, designed to ensure that, on average, only a certain proportion (e.g., 5%) of the reported significant results are false positives [90] [88].

The Critical Role of Dependency in Splicing Data

A fundamental assumption underlying many multiple testing corrections, including the Bonferroni method, is the independence of tests [85]. This assumption is profoundly violated in biological data, especially in alternative splicing analysis.

• Biological Dependency: Splicing events within a gene are not independent. The inclusion or exclusion of one exon can influence the splicing patterns of neighboring exons through mechanisms like exon definition and splice site competition [2]. Furthermore, genes operate in coordinated pathways, leading to correlated expression and splicing changes across gene families [85].
• Consequences for FDR: While the BH procedure is mathematically proven to control FDR under positive dependency, this guarantee pertains to the expected value of the False Discovery Proportion (FDP) across many hypothetical experiments [90]. In any single experiment with strongly correlated features—such as a set of co-regulated exons—the actual FDP can be highly variable and much larger than the nominal FDR level [90]. Recent studies have shown that in datasets with strong dependencies (e.g., metabolomics, correlated methylation sites), standard BH correction can, in a subset of experiments, yield a surprisingly high number of false positives, even when all null hypotheses are true [90].
• Empirical Evidence: This phenomenon has been demonstrated in RNA-seq data. When group labels were shuffled in a differential expression analysis (creating a true null scenario), the application of BH correction still resulted in some experiments where a high proportion of genes were falsely called significant. These false-positive genes were found to be highly correlated with each other [90]. This directly translates to splicing analysis, where dependencies are even more pronounced.

Experimental and Computational Protocols for Alternative Splicing Analysis with FDR Control

Protocol 1: Differential Splicing Analysis with rMATS and BH Correction

This protocol is adapted from a standard bulk RNA-seq splicing workflow [62].

I. RNA Sequencing and Primary Data Generation

• RNA Extraction & Library Prep: Isolate total RNA using a validated kit (e.g., RNeasy Plus). Prepare stranded mRNA sequencing libraries using a commercial kit (e.g., KAPA mRNA HyperPrep). Include a minimum of four biological replicates per condition to ensure statistical power [62].
• Sequencing: Sequence libraries on an Illumina platform (e.g., NovaSeq) using paired-end 150 bp reads. Aim for a minimum depth of 30-50 million reads per sample to adequately capture junction-spanning reads.
• Alignment: Align trimmed FASTQ reads to a reference genome (e.g., GRCh38) using a splice-aware aligner such as STAR. This step is critical for accurate identification of splice junctions.

II. Splicing Quantification with rMATS

• Input Preparation: Provide rMATS with the coordinate-sorted BAM files for all replicates, grouped by condition.
• Event Identification: Run rMATS (v4.1.0 or higher) to identify and quantify five basic types of alternative splicing events: Skipped Exon (SE), Alternative 5' Splice Site (A5SS), Alternative 3' Splice Site (A3SS), Mutually Exclusive Exons (MXE), and Retained Intron (RI).
• Output Metrics: For each event, rMATS calculates:
  • PSI (Ψ): Percent Spliced-In, the ratio of reads supporting the inclusion isoform.
  • ΔPSI: The difference in PSI between two experimental conditions.
  • P-value: The raw probability from a likelihood ratio test that ΔPSI ≠ 0.
  • FDR: The Benjamini-Hochberg adjusted q-value for the event, correcting for all tests performed within that event type [62].

III. Result Interpretation and Filtering

• Thresholding: Apply dual filters to identify high-confidence differential splicing events (DSEs). Common thresholds are |ΔPSI| > 0.1 (or 10%) and FDR < 0.05.
• Visualization: Inspect significant events using Integrated Genome Viewers or generate sashimi plots to validate read coverage and junction support.

Protocol 2: Handling Heterogeneous Data with MAJIQ HET

For large, heterogeneous datasets (e.g., from public repositories like GTEx with many donors and tissues), traditional group-wise models can fail. The MAJIQ v2 framework introduces a heterogeneous (HET) test statistic to address this [15].

I. Splice Graph Construction and Quantification

• Build Splice Graph: Run the MAJIQ builder to create a unified splice graph for each gene from all samples. This graph incorporates both annotated and de novo (unannotated) splice junctions, which are common in disease states [15].
• Define LSVs: Identify Local Splicing Variations (LSVs), which represent splits in the splice graph entering or leaving a reference exon.
• Quantify Heterogeneously: Use the MAJIQ HET quantifier. Instead of assuming a shared PSI per condition, it calculates a posterior distribution for PSI (𝝍) for each sample individually [15].

II. Non-Parametric Differential Splicing Test

• Apply Rank-Based Test: For each LSV, compare the distributions of per-sample PSI estimates between two conditions using a robust, non-parametric test (e.g., Mann-Whitney U or a custom TNOM test) [15].
• Calculate Confidence: Compute the probability that the absolute ΔPSI exceeds a user-defined minimum threshold (e.g., P(|ΔΨ| > 0.1)) [15].
• FDR Control: Apply the Benjamini-Hochberg procedure to the p-values (or equivalent scores) generated from the non-parametric tests across all LSVs to control the FDR.

Table 2: Comparison of Splicing Analysis Tools & FDR Considerations

Tool	Core Methodology	Unit of Analysis	Handles Dependencies	FDR Control Method	Best For
rMATS [62]	Likelihood ratio test on junction counts.	Pre-defined AS event types (SE, MXE, etc.).	No explicit modeling.	Standard Benjamini-Hochberg.	Standard designs with clear biological replicates.
MAJIQ v2 [15]	Bayesian modeling of splice graphs, non-parametric HET test.	Local Splicing Variations (LSVs).	Yes, via splice graphs and sample-specific PSI.	Benjamini-Hochberg on HET test outputs.	Large, heterogeneous datasets (e.g., GTEx).
LeafCutter	Clustering of intron excision counts.	Splicing clusters (introns).	Yes, via intron clustering.	Permutation-based FDR.	Discovery of novel and complex splicing events.

Table 3: Research Reagent Solutions for Bulk RNA-seq Splicing Analysis

Item / Resource	Function / Purpose	Example / Specification
RNA Isolation Kit	To extract high-integrity total RNA, free of genomic DNA, crucial for accurate transcriptome representation.	RNeasy Plus Micro/Mini Kit (Qiagen) [62].
Stranded mRNA Library Prep Kit	To generate sequencing libraries that preserve strand-of-origin information, essential for distinguishing overlapping transcripts.	KAPA mRNA HyperPrep Kit [62].
Splice-Aware Aligner	To accurately map sequencing reads across exon-exon junctions, the foundation of all splicing analysis.	STAR (Spliced Transcripts Alignment to a Reference).
Splicing Quantification Software	To identify splicing events, calculate inclusion levels (PSI), and perform statistical testing for differential splicing.	rMATS [62], MAJIQ [15], LeafCutter.
FDR Control Package	To implement multiple testing corrections on p-values generated by analysis tools.	p.adjust function in R (method="BH"), statsmodels.stats.multitest.fdrcorrection in Python.
Reference Transcriptome	A comprehensive set of annotated transcripts and splice junctions for alignment and quantification.	GENCODE, Ensembl (e.g., GRCh38 release) [62].

Advanced Considerations and Best Practices

Given the dependency structure in splicing data, blind reliance on a nominal FDR threshold (e.g., 5%) can be misleading [90]. Adopting a layered, skeptical approach is recommended.

• Employ Negative Controls: Where possible, incorporate synthetic null data or negative control samples (e.g., samples from the same condition with shuffled labels) into the analysis pipeline. This can help empirically estimate the false discovery proportion specific to your dataset's correlation structure [90].
• Prioritize Robust Events: Place greater confidence in splicing events that are supported by large ΔPSI magnitudes and high read coverage. Events with a ΔPSI of 0.4 are less likely to be false positives than those at the 0.1 cutoff, even with identical FDR values.
• Independent Validation: Plan for orthogonal validation of key findings using techniques such as RT-PCR with capillary electrophoresis or nanopore long-read sequencing. This remains the gold standard for confirming splicing changes.
• Contextualize FDR Outputs: Understand that an FDR of 5% does not guarantee that exactly 5% of your reported DSEs are false; rather, it is a statement about the long-run average. In a given experiment with correlated events, the actual proportion may differ [90]. Reporting both FDR and measures of effect size (ΔPSI) is essential.
• Explore Dependency-Aware Methods: For mission-critical analyses, consider methods that explicitly account for correlation, such as permutation-based FDR estimation (used in eQTL studies to handle linkage disequilibrium) [90] or the heterogeneous testing in MAJIQ v2 [15].

Dealing with Low-Abundance Isoforms and Complex Splicing Variations

Within the broader thesis on alternative splicing analysis in bulk RNA sequencing research, the precise identification and quantification of low-abundance isoforms represent a critical frontier. These rare transcripts, often constituting less than 10% of a gene's total expression, can encode functionally distinct proteins or regulatory non-coding RNAs pivotal to cellular differentiation, disease progression (e.g., cancer and neurodegeneration), and drug response [91]. However, their study is confounded by technical limitations: the short reads from dominant sequencing platforms like Illumina struggle to uniquely map to repetitive or homologous regions, leading to ambiguous mapping and quantification inaccuracy [91] [92]. Furthermore, the dynamic range of expression in a bulk sample masks cell-type-specific splicing events, conflating signals from multiple cell populations [91] [93].

The strategic framework for addressing these challenges rests on three pillars: 1) Enhanced Experimental Design, prioritizing sequencing depth (>100M reads per sample for complex genomes) and long-read or strand-specific protocols to capture full-length transcript context [91] [92]; 2) Advanced Computational Quantification, employing probabilistic models that resolve read multi-mapping to estimate isoform abundance accurately [94] [95]; and 3) Multi-Omics Integration, correlating splicing variations with proteomic data or functional screens to validate biological significance [96] [97]. This integrated approach is directly aligned with national research priorities, such as the "Empowering Drug Innovation through RNA Basic Research" major project, which emphasizes the discovery of novel RNA diagnostic markers and therapeutic targets through multi-modal data analysis [96].

Methodological Comparison: Tools for Isoform Detection and Quantification

Selecting the appropriate bioinformatics tool is paramount for reliable isoform analysis. The tools differ fundamentally in their algorithms, inputs, and optimal use cases.

Table 1: Comparison of Key Tools for Transcript-Level Analysis from RNA-Seq Data

Tool	Core Algorithm	Input Requirements	Key Advantages	Limitations / Best For	Typical Use Case in Splicing Analysis
Salmon [94]	Alignment-free, quasi-mapping using k-mers & a rich statistical model.	Raw FASTQ files + transcriptome FASTA (decoy-aware index).	Extremely fast, accounts for GC/sequence bias, accurate for lowly-expressed transcripts.	Does not generate genomic BAM files; quantification relies on a provided transcriptome.	Rapid, large-scale quantification of known and novel transcripts from large cohorts.
RSEM (RNA-Seq by Expectation-Maximization) [95]	Alignment-based EM algorithm for estimating maximum likelihood abundances.	BAM file aligned to transcriptome or genome + GTF annotation.	Very accurate, integrates well with aligners (Bowtie2, STAR), generates expected read counts and credibility intervals.	Slower than alignment-free methods; dependent on the accuracy of initial alignment.	High-precision quantification where alignment information and visual validation (BAM files) are required.
StringTie [98]	Network flow algorithm for de novo transcript assembly and quantification.	Genome-aligned BAM file (e.g., from HISAT2/STAR) + reference GTF (optional).	Can discover novel isoforms not in the annotation; performs assembly and quantification simultaneously.	Assembly can be computationally intensive; novel isoforms require careful validation.	Discovery-focused studies aiming to identify previously unannotated splicing variants and fusion transcripts.
SCASL (Single-Cell Alternative Splicing Landscape) [93]	k-NN imputation and spectral clustering specifically for sparse single-cell data.	Single-cell RNA-seq count matrices (e.g., from CellRanger).	Resolves splicing heterogeneity at single-cell resolution; defines cell states based on splicing profiles.	Exclusively for single-cell data; requires moderate to high sequencing depth per cell.	Defining tumor subclones, developmental trajectories, or immune cell states based on cell-specific splicing patterns.

Detailed Experimental Protocol: From Sample to Isoform-Calling

This protocol outlines a robust, tiered strategy for identifying and validating low-abundance splicing variants from bulk tissue.

Phase 1: Sample Preparation & High-Depth Sequencing

• RNA Extraction & QC: Isolate total RNA using a method that preserves integrity (RIN > 8). Treat samples with DNase I. Use a ribosomal RNA depletion kit (not poly-A selection) to retain both polyadenylated and non-polyadenylated transcripts, including some circular RNAs [97].
• Strand-Specific Library Construction: Prepare libraries using a stranded mRNA-seq protocol [92]. This preserves the information of the transcript of origin, crucial for accurately determining the expression of overlapping genes on opposite strands and for correct isoform assembly.
• Sequencing Platform Selection: For discovery, use an Illumina short-read platform (e.g., NovaSeq) to generate a minimum of 150 million paired-end (PE) 150bp reads per sample. This depth is necessary to capture reads spanning low-abundance splice junctions [91]. For validation and to resolve complex loci, supplement with long-read sequencing (e.g., PacBio Revio or Oxford Nanopore) on a subset of samples to obtain full-length isoform sequences [91].

Phase 2: Bioinformatic Processing & Quantification

• Quality Control & Trimming: Use FastQC for initial assessment. Trim adapters and low-quality bases with Trimmomatic (parameters: SLIDINGWINDOW:4:20 MINLEN:25) [98].
• Reference-Based Alignment & Assembly: Align reads to the reference genome using a splice-aware aligner like STAR or HISAT2 [98]. Use the aligned BAM file as input for StringTie2 in a reference-guided assembly mode. This step will generate a merged transcriptome annotation (GTF file) that includes both known and novel isoforms [98].
• Transcript Quantification: Build a Salmon index using the merged transcriptome FASTA (generated from the GTF). Quantify each sample against this index using Salmon in alignment-free mode with bias correction flags (--seqBias --gcBias). This combines the discovery power of assembly with the speed and accuracy of lightweight quantification [94].
• Differential Splicing Analysis: Import Salmon's transcript-level counts into R/Bioconductor using the tximport package. Use specialized differential splicing/testing tools like DEXSeq (for differential exon usage) or spliceVar to identify isoforms or splicing events significantly altered between conditions.

Phase 3: Validation & Functional Prioritization

• PCR Validation: Design primers spanning novel or low-abundance splice junctions identified. Perform RT-PCR followed by Sanger sequencing or use digital droplet PCR (ddPCR) for absolute, sensitive quantification.
• Functional Correlation: Integrate isoform expression data with proteomics data (if available) to check for correlation. Perform gene ontology (GO) and pathway enrichment analysis on genes producing significant low-abundance isoforms to infer potential functional impact.
• Single-Cell Contextualization: For tissues, leverage public or generate new single-cell RNA-seq data. Use tools like SCASL to investigate whether the low-abundance isoform is specific to a rare but functionally important cell subpopulation, which would explain its low signal in bulk data [93].

Visualizing Experimental Strategy and Tool Selection

Low-Abundance Isoform Analysis Workflow

Tool Selection Logic for Splicing Analysis

Table 2: Research Reagent Solutions for Low-Abundance Isoform Analysis

Category	Item / Kit	Function / Purpose	Key Considerations
Sample Prep	Ribosomal RNA Depletion Kits (e.g., NEBNext rRNA Depletion)	Removes abundant rRNA, enriching for mRNA, lncRNA, and other transcripts, providing a more uniform coverage crucial for detecting low-abundance species.	Preferred over poly-A selection for broader transcriptome coverage, especially for non-polyadenylated or partially degraded RNAs [97] [92].
	RNase Inhibitors & High-RIN Extraction Kits	Preserves RNA integrity during isolation, minimizing degradation that can obscure or destroy rare transcript signals.	Essential for working with sensitive tissues or low-input samples. Aim for RNA Integrity Number (RIN) > 8.
Library Prep	Stranded mRNA Library Prep Kits (e.g., Illumina Stranded mRNA Prep)	Preserves strand-of-origin information during cDNA synthesis, allowing unambiguous determination of overlapping antisense transcripts and correct isoform assembly [92].	A standard for modern splicing studies.
	Long-Read cDNA Prep Kits (e.g., PacBio Iso-Seq or ONT cDNA kits)	Generates full-length cDNA for long-read sequencing, resolving complex splice variants and novel isoforms unambiguously [91].	Used for targeted validation and discovery on a subset of samples due to higher cost and lower throughput.
Sequencing	Illumina NovaSeq S4 Flow Cell	Provides the ultra-high sequencing depth (>100M PE reads/sample) required to achieve sufficient coverage for statistically robust detection of low-abundance splice junctions [91].	The workhorse for deep discovery phases. Cost-effective per gigabase.
	PacBio Revio or Sequel IIe System	Enables HiFi long-read sequencing for definitive isoform resolution, phasing of multiple variants, and discovery of complex structural variations [91].	Used as a complementary technology to short-read data.
Validation	Digital Droplet PCR (ddPCR) Supermix	Provides absolute, sensitive quantification of specific splice junctions or isoforms without the need for a standard curve, ideal for validating low-abundance targets.	Higher sensitivity and precision than qPCR for rare targets.
	High-Fidelity DNA Polymerase for PCR (e.g., PrimeSTAR GXL)	Amplifies specific, potentially long cDNA fragments spanning splice junctions for Sanger sequencing validation of novel isoforms.	High fidelity reduces PCR errors, crucial for confirming sequence.

Incorporating De Novo and Unannotated Splicing Events

Within the broader thesis of alternative splicing (AS) analysis in bulk RNA sequencing research, a critical frontier is the move beyond canonical, annotated transcriptomes. Relying solely on existing annotations restricts discovery, as a significant portion of biologically relevant and condition-specific splicing variation remains uncharacterized. De novo and unannotated splicing events—junctions, exons, and intron retentions not present in reference databases—represent a substantial layer of transcriptomic complexity. Their systematic incorporation is not merely a technical improvement but a conceptual necessity for a complete understanding of gene regulation in development, cell-type specification, and disease pathogenesis. This article provides detailed application notes and protocols for detecting, quantifying, and validating these novel events, framing them as essential components of a modern, comprehensive alternative splicing analysis workflow.

Quantitative Evidence: The Prevalence and Impact of Unannotated Splicing

The biological imperative to incorporate unannotated events is grounded in their demonstrated prevalence and functional impact. Large-scale analyses of public RNA-seq data reveal that a stable, significant fraction of expressed splice junctions are absent from major annotations.

Table 1: Prevalence of Unannotated Splicing Junctions in Human Transcriptome

Metric	Finding	Data Source & Scope	Implication for Research
Junctions in ≥1000 samples	18.6% (56,861 junctions) were not annotated [99].	21,504 human RNA-seq samples from SRA [99].	Canonical annotations miss a large, recurrent set of splicing variations.
Junctions in ≥40% of samples	>99.8% were annotated [99].	Same as above [99].	Highly common junctions are well-documented; novel events are often context-specific.
Increase in detected differential splicing	Accounting for complex/unannotated events yielded a >30% increase in detected differentially spliced events [15].	Comparison across mouse tissues [15].	Incorporating novel events dramatically increases analytical sensitivity and discovery.
Cell-type-specific novel events	In neuronal subclasses, up to 30% of detected differential splicing events were novel [100].	Analysis of mouse brain single-cell datasets [100].	Novel splicing is particularly prevalent in complex tissues and may define specialized cell states.

Incorporating these events is especially critical in disease contexts like cancer and neurodegeneration, where aberrant splicing is a hallmark, and in defining the molecular identity of cell types, where splicing can be a more precise marker than gene expression alone [15] [101].

Core Computational Workflow and Protocols

The following protocol outlines a robust, annotation-aware yet annotation-free capable pipeline for bulk RNA-seq AS analysis, based on state-of-the-art tools like MAJIQ v2 and principles from LeafCutter/scQuint [15] [100].

Protocol 1: De Novo Splicing Analysis from Bulk RNA-Seq Alignment

Objective: To identify and quantify both annotated and unannotated local splicing variations (LSVs) from raw RNA-seq data.

Input: Paired-end FASTQ files for all samples; reference genome (FASTA); optional gene annotation (GTF/GFF3). Output: Quantified LSVs, including de novo junctions, with Percent Spliced-In (Ψ, PSI) values and differential splicing probabilities.

Step-by-Step Procedure:

• Splice-Aware Alignment:

  • Align reads to the reference genome using a splice-aware aligner (e.g., STAR [18]).
  • Critical Parameter: Use --twopassMode Basic to discover novel junctions. Set --limitSjdbInsertNsj sufficiently high (e.g., 2000000) to accommodate novel junction discovery in large genomes.
  • Output coordinate-sorted BAM files for each sample.

• Build or Update the Splice Graph:

  • Use MAJIQ Builder to process BAMs and create a unified splice graph per gene [15].
  • Provide optional annotation (GTF) as a guide but do not restrict to it.
  • The algorithm identifies all possible junctions and intron retentions from read coverage, adding de novo elements (unannotated junctions, exons) to the graph [15].
  • Filtering: Apply user-defined filters (e.g., --min-reads per junction, --min-samples per group) to remove low-confidence, likely artifacial events [15].
  • This step is incremental; new samples can be added without reprocessing the entire dataset [15].

• Quantify Local Splicing Variations (LSVs):

  • Run MAJIQ Quantifier on the built graph [15].
  • It defines LSVs as splits in the splice graph entering or leaving a reference exon, capturing complex events beyond simple cassette exons [15].
  • For each LSV junction in each sample, it calculates a posterior distribution for its PSI (Ψ) using a Bayesian model that accounts for read counts and distribution [15].
  • For heterogeneous datasets or large sample groups, use the MAJIQ HET module, which quantifies PSI per sample and uses robust rank-based statistics (e.g., Mann-Whitney U) for group comparisons [15].

• Differential Splicing Analysis:

  • Specify sample groups (e.g., Control vs. Treated).
  • MAJIQ computes the posterior for the change in PSI (ΔΨ) between conditions.
  • The primary output is the probability that |ΔΨ| > a user-defined threshold (e.g., C=0.2), denoted as P(|ΔΨ| > C). Events with high probability (e.g., >0.95) are high-confidence differential splicing events [15].

• Visualization and Interpretation (VOILA v2):

  • Use VOILA v2 to visualize complex splicing graphs, PSI distributions across hundreds of samples, and ΔΨ values for significant LSVs [15].
  • The Modulizer algorithm can optionally classify LSVs into coherent AS modules (e.g., cassette exon, alternative 3'/5' site, intron retention) to simplify interpretation [15].

Diagram 1: Computational Workflow for De Novo Splicing Analysis. The pipeline progresses from raw data to high-confidence differential splicing events involving unannotated junctions, integrating key steps for splice graph construction and robust statistical testing.

Experimental Validation of Novel Splicing Events

Computational predictions of de novo splicing require empirical validation. PCR-based methods remain the gold standard for gene-specific confirmation.

Protocol 2: RT-PCR and Capillary Electrophoresis Validation

Objective: To experimentally confirm the existence and quantify the relative abundance of predicted novel splice isoforms.

Step-by-Step Procedure:

• Primer Design:

  • Design primers in exons flanking the predicted novel event (e.g., spanning an unannotated junction or retained intron).
  • Ensure one primer is fluorescently labeled (e.g., 6-FAM) for capillary detection.
  • For novel junctions, design one primer overlapping the predicted exon-exon boundary to ensure specificity.

• cDNA Synthesis and PCR:

  • Use high-quality, DNase-treated total RNA.
  • Perform reverse transcription with random hexamers and a high-fidelity reverse transcriptase.
  • Run PCR with a high-fidelity polymerase, using a low cycle number (e.g., 25-30) to remain in the exponential phase and allow semi-quantitative comparison [1].
  • Include a no-template control and a positive control (sample where the event was computationally predicted).

• Product Separation and Detection:

  • Option A - Capillary Electrophoresis (Recommended): Mix PCR products with size standard and run on a genetic analyzer (e.g., ABI PRISM). This method provides single-base-pair resolution, precise sizing, and quantitative fluorescence peak areas for each isoform [1].
  • Option B - Agarose Gel Electrophoresis: Separate PCR products on a high-percentage agarose gel (2-3%). Extract bands, purify, and submit for Sanger sequencing to definitively confirm the novel junction sequence [1].

• Quantitative Analysis:

  • For capillary data, the relative abundance of an isoform can be estimated as the percentage of its peak area relative to the total peak area for all isoforms from that gene in the sample [1].
  • Compare these percentages between experimental conditions to validate predicted ΔΨ trends.

Diagram 2: Experimental Validation Pipeline for De Novo Splicing. Two primary detection paths lead from PCR product to confirmed novel isoform, with capillary electrophoresis enabling direct quantification.

Comparative Analysis of Splicing Detection Methods

Selecting the appropriate computational tool is crucial. The table below compares key methods that handle unannotated events, adapted for bulk RNA-seq analysis.

Table 2: Comparison of Computational Methods for De Novo Splicing Detection

Method	Core Approach	Handles Unannotated Events?	Strengths for De Novo Analysis	Key Considerations
MAJIQ v2 [15]	Local Splicing Variations (LSVs) via splice graph.	Yes. Explicitly adds de novo junctions/exons to graph.	Incremental build; HET module for heterogeneous data; excellent visualization (VOILA2).	Can be computationally intensive for very large graphs.
LeafCutter / scQuint [100]	Clustering of alternative intron usage.	Yes. Annotation-free; detects events from intron clusters.	Highly robust to coverage biases; scalable.	Focuses on intron-centric events; may miss exon-centered complexity.
rMATS [101]	Statistical modeling of junction counts for defined event types.	Limited. Primarily quantifies pre-defined, annotated event types.	Fast, standardized for common event types.	Not designed for discovery of unannotated event structures.
Whippet / MicroExonator [101]	Splicing node quantification (e.g., core exons, donors).	Can detect unannotated nodes from reference.	Efficient for large datasets; good for micro-exons.	Relies on a reference for node definition.

The Scientist's Toolkit: Essential Research Reagents and Materials

Table 3: Key Reagent Solutions for Splicing Discovery and Validation

Reagent / Material	Function in Protocol	Specific Application Note
High-Fidelity Reverse Transcriptase	Converts RNA to cDNA with minimal enzyme-introduced errors, preserving true splice variant sequences.	Essential for all downstream PCR validation steps. Use with random hexamers for unbiased representation [1].
Hot-Start High-Fidelity DNA Polymerase	Amplifies specific cDNA targets with high accuracy and yield, minimizing nonspecific products during PCR validation.	Critical for amplifying low-abundance novel isoforms. The hot-start feature improves specificity [1].
Fluorophore-labeled PCR Primers	Enables detection and precise quantification of PCR amplicons via capillary electrophoresis.	Typically 5'-labeled with 6-FAM, HEX, or similar. Allows multiplexing and high-resolution fragment analysis [1].
Size Standard for Capillary Electrophoresis	Provides an internal ladder for accurate sizing of PCR amplicons to the exact base pair.	Crucial for distinguishing isoforms with small size differences (e.g., alternative 3' splice sites differing by few bp) [1].
DNase I (RNase-free)	Removes contaminating genomic DNA from RNA preparations to prevent false-positive signals from unspliced introns.	A mandatory step before cDNA synthesis, especially when primers are located in different exons [1].
Splice-Aware Aligner Software (STAR)	Aligns RNA-seq reads to the genome, correctly spanning introns and identifying canonical and novel splice junctions.	The --twopassMode parameter is key for de novo junction discovery [18].
Reference Splicing Annotation (GTF)	Provides a baseline of known transcripts to guide but not limit splice graph construction.	Using a comprehensive set (e.g., GENCODE) improves initial graph building but should not be used as a filter [15] [99].

Diagram 3: Splice Graph Incorporating De Novo Elements. A visual representation of how novel exons, junctions, and splice site variants (dashed lines) are integrated with the annotated transcript structure (solid lines) to form a complete splicing graph.

Best Practices for RNA Extraction, Library Preparation, and Sequencing Depth

This document outlines best practices and detailed protocols for bulk RNA sequencing (RNA-seq) experiments, specifically tailored for research focused on alternative splicing analysis. Accurate detection and quantification of splicing events, which are critical for understanding gene regulation in development and disease, impose stringent requirements on every step of the workflow—from RNA integrity and library construction to sequencing depth [102]. The recommendations herein are designed to ensure the generation of high-quality, reproducible data capable of revealing subtle isoform-level changes, thereby supporting robust conclusions within a broader thesis on splicing dynamics.

RNA Extraction and Quality Control

The fidelity of any RNA-seq experiment is fundamentally limited by the quality of the input RNA. For splicing analysis, where the integrity of full-length transcripts is paramount, rigorous extraction and quality assessment are non-negotiable.

Core Principles and Best Practices

• RNase Inhibition: RNases are ubiquitous and stable. Use an RNase-free workspace decontaminated with specialized solutions (e.g., RNaseZap), and employ RNase-free consumables [103] [104].
• Immediate Stabilization: Upon collection, stabilize RNA immediately to halt degradation. Methods include flash-freezing in liquid nitrogen, homogenization in chaotropic lysis buffers (e.g., TRIzol), or immersion in stabilization reagents like RNAlater [103].
• Appropriate Isolation Method: Choose a kit matched to your sample type.
  • Fresh cells/tissues: Silica-membrane column kits (e.g., PureLink RNA Mini Kit) offer a good balance of ease, quality, and yield [103].
  • Low-input samples (<10,000 cells): Use dedicated low-input kits (e.g., PicoPure RNA Isolation Kit) [105] [19].
  • Fibrous or lipid-rich tissues: Acid-guanidinium-phenol-chloroform (e.g., TRIzol) extraction may be more effective [103].
  • Formalin-Fixed Paraffin-Embedded (FFPE) tissue: Use kits specifically designed to reverse cross-links, providing higher yields of amplifiable RNA [104].
• DNase Treatment: Perform on-column DNase digestion to remove genomic DNA contamination, which is crucial for accurate quantification and for applications using intron-spanning primers [103].
• Proper Storage: Aliquot RNA and store at -80°C for long-term preservation. Avoid repeated freeze-thaw cycles [103] [104].

Quality Assessment Protocol

Do not proceed to library preparation without verifying RNA quality and quantity using the following metrics [103] [104] [106].

Quantification and Purity:

• Use a UV-Vis spectrophotometer (e.g., NanoDrop). Acceptable purity ratios are:
  • A260/A280: ~1.8–2.0 (indicates low protein contamination).
  • A260/A230: >2.0 (indicates low organic compound contamination).
• For more accurate quantification, especially with low-concentration samples, use fluorescence-based assays (e.g., Qubit RNA HS Assay).

Integrity Assessment:

• Analyze RNA integrity using capillary electrophoresis (e.g., Agilent Bioanalyzer or TapeStation).
• Key metrics:
  • RNA Integrity Number (RIN): A score ≥ 7.0 is generally recommended for standard poly(A)-selection protocols [19] [106]. Splicing analysis benefits from even higher RIN (>8.0).
  • DV200: The percentage of RNA fragments > 200 nucleotides. This metric is particularly valuable for degraded samples (e.g., FFPE). A DV200 > 70% is ideal, while protocols can be adapted for samples with DV200 between 30–70% [107].

Table 1: RNA Quality Metrics and Their Implications for Library Prep

Metric	Ideal Value	Minimum Value	Implication for Splicing Analysis
RIN	8.0 - 10.0	≥ 7.0 [106]	High integrity ensures full-length transcripts, enabling accurate junction mapping.
DV200	> 70%	> 50% [107]	Critical for degraded samples; low DV200 requires more reads to cover spliced transcripts.
A260/A280	1.9 - 2.1	≥ 1.8	Low ratio suggests protein/phenol contamination, which can inhibit downstream enzymatic steps.
Total RNA Input	100 - 1000 ng	10 - 300 ng [105]	Low input demands specialized kits and may increase duplicate rates, affecting complexity.

Library Preparation for Splicing Analysis

Library preparation determines what RNA species are sequenced and how accurately their original orientation and structure are preserved. For splicing studies, stranded, rRNA-depleted libraries are often optimal.

Key Strategic Decisions

• Poly(A) Selection vs. Ribosomal RNA (rRNA) Depletion: Poly(A) selection captures mRNA with intact poly-adenylated tails, but it biases against non-polyadenylated transcripts and is unsuitable for degraded RNA (low RIN/DV200). rRNA depletion removes ribosomal sequences (≈80% of total RNA) and preserves other RNA biotypes, making it superior for samples with lower integrity or for studying non-coding RNAs [106]. For example, the Polaris Depletion kit targets both rRNA and globin mRNA [105].
• Strandedness: Always use a stranded library preparation protocol. This preserves the information about which DNA strand was transcribed, allowing you to accurately distinguish overlapping genes on opposite strands and to determine the correct orientation of splicing junctions, which is essential for isoform-level analysis [106].
• Unique Molecular Identifiers (UMIs): For low-input or highly multiplexed experiments where PCR duplication rates are high, incorporate UMIs. They allow bioinformatic collapse of reads originating from the same original molecule, improving quantification accuracy [108] [107].

Detailed Library Preparation Protocol

The following protocol is adapted from the JAX_DPC method using the Watchmaker RNA Library Prep Kit with Polaris Depletion, suitable for bulk RNA-seq from varied sample types like organoids and iPSCs [105].

Materials: Watchmaker RNA Library Prep Kit, Polaris Depletion Kit, IDT for Illumina TruSeq DNA UD Indexes, RNase-free reagents, magnetic stand, thermal cycler.

Procedure:

• Input RNA: Begin with high-quality total RNA. A typical input is 300 ng, but the protocol can accommodate a range [105].
• rRNA/Globin Depletion (Polaris Depletion):
  • Incubate the RNA with depletion probes that hybridize to rRNA and globin mRNA sequences.
  • Add RNase H to digest the RNA in the RNA-DNA hybrid complexes.
  • Purify the remaining, enriched RNA using magnetic beads.
• First-Strand cDNA Synthesis: Using random primers, reverse transcribe the depleted RNA into first-strand cDNA.
• Second-Strand Synthesis & End Repair: Synthesize the second strand to create double-stranded cDNA. Repair ends to create blunt fragments.
• Adapter Ligation: Ligate indexed sequencing adapters to the cDNA fragments. Using dual-indexed adapters (e.g., TruSeq DNA UD Indexes) is critical for sample multiplexing and demultiplexing [105].
• Library Amplification: Perform a limited-cycle PCR to enrich for adapter-ligated fragments and add full adapter sequences required for sequencing.
• Library Quality Control and Quantification:
  • Assess library fragment size distribution using a Bioanalyzer (expect a broad peak centered around 300-500 bp).
  • Quantify libraries accurately by qPCR (e.g., Kapa Library Quantification Kit) for precise pooling and loading onto the sequencer.

Diagram 1: RNA-seq Library Preparation Workflow

Sequencing Depth and Experimental Design

Sequencing depth is the single most important factor determining the power to detect differential splicing events, especially for low-abundance transcripts [102].

Determining Optimal Depth for Splicing

The standard depth for differential gene expression (25-40 million reads) is insufficient for comprehensive splicing analysis. Splicing detection requires coverage across exon-exon junctions, which are less abundant than exonic reads.

• Foundation Data: A 2025 study on Mendelian disorders systematically evaluated depth, finding that pathogenic splicing abnormalities undetectable at 50 million reads became clearly apparent at 200 million reads and were fully characterized at 1 billion reads [102].
• Saturation Point: While gene detection saturates around 1 billion reads, isoform-level detection continues to improve with increasing depth, indicating a higher complexity [102].
• General Guideline: For human splicing-focused studies, aim for a minimum of 100 million paired-end reads per sample. Complex discoveries or low-abundance event detection may require 200 million reads or more [107].

Table 2: Recommended Sequencing Depth by Research Goal

Research Goal	Recommended Minimum Depth (Paired-end Reads)	Read Length	Justification
Differential Gene Expression	25 - 40 million [107]	2x75 bp	Cost-effective, stabilizes fold-change estimates for medium-high abundance genes.
Alternative Splicing / Isoform Detection	100 - 200 million [102] [107]	2x100 bp or 2x150 bp	Required to achieve sufficient coverage across splice junctions for accurate quantification of isoforms.
Fusion Gene Detection	60 - 100 million [107]	2x100 bp	Needs enough reads spanning breakpoints to identify chimeric junctions with confidence.
Variant Calling / Allele-Specific Expression	≥ 100 million [107]	2x100 bp	High depth reduces sampling error for accurate allele frequency estimation.

Integrating Depth with RNA Quality and Experimental Design

• Degraded RNA: If RNA integrity is suboptimal (DV200 30-50%), increase sequencing depth by 25-50% to compensate for reduced library complexity and increased duplication [107].
• Read Length: Longer reads (2x100 bp or 2x150 bp) improve the mappability of reads spanning splice junctions and are strongly recommended over 2x75 bp for splicing analysis [107].
• Replicates: Biological replicates are non-negotiable for statistical power. For splicing, at least 3-5 biological replicates per condition are recommended to account for biological variability in isoform usage.
• Controlling for Confounds: Document all metadata (sample origin, extraction batch, library prep date, sequencer lane) meticulously. Use this metadata during statistical analysis to control for batch effects, which are a major source of false positives [109].

Diagram 2: Experimental Design Logic for Splicing Studies

The Scientist's Toolkit: Essential Reagents & Materials

Table 3: Key Research Reagent Solutions for Bulk RNA-seq

Reagent / Kit	Primary Function	Key Consideration for Splicing Analysis
TRIzol Reagent	Guanidinium-phenol-based total RNA isolation.	Robust for difficult tissues (high fat, RNase); handles varied sample types [103].
PureLink RNA Mini Kit	Silica-membrane column purification of total RNA.	Easy, reliable method for most fresh samples; includes on-column DNase option [103].
PicoPure RNA Isolation Kit	RNA purification from low cell counts (<10,000 cells).	Essential for limited samples like sorted cell populations [105] [19].
RNase-free DNase Set	Removal of contaminating genomic DNA.	Critical step to prevent false positives, especially for genes with pseudogenes [103].
Watchmaker RNA Lib Prep Kit w/ Polaris Depletion	rRNA/globin depleted, stranded library prep.	Preserves full transcriptome information; ideal for medium-quality RNA and splicing [105].
NEBNext Ultra II Directional RNA Library Prep	Poly(A)-selected, stranded library prep.	Standard for high-quality mRNA sequencing; good for high-RIN samples [19].
IDT for Illumina TruSeq DNA UD Indexes	Dual-index adapters for sample multiplexing.	Allows high-level multiplexing, reducing per-sample cost and batch effects [105].
RNaseZap / RNaseX	Surface decontamination solution.	Foundational step to protect RNA integrity during all handling steps [103] [104].

Validating Splicing Events and Benchmarking Computational Methods

Within the broader scope of a thesis focused on alternative splicing (AS) analysis via bulk RNA sequencing (RNA-seq), the discovery and quantification of splicing events represent a critical computational challenge [2]. High-throughput RNA-seq enables transcriptome-wide profiling but introduces inherent uncertainties in isoform-level quantification, particularly for low-abundance transcripts or complex splicing patterns [19]. Computational tools like rMATS can calculate metrics such as the Percent Spliced-In (PSI) to identify differential splicing events between conditions [62]. However, these in silico predictions require rigorous, orthogonal experimental validation to confirm their biological reality and accuracy before downstream mechanistic or clinical conclusions can be drawn.

This document establishes multiplex Reverse Transcription-Polymerase Chain Reaction (RT-PCR) coupled with capillary electrophoresis (CE) as a gold standard methodology for this essential validation step. The integration of RT-PCR’s high sensitivity and specificity with CE’s superior fragment resolution and multiplexing capability provides a robust, quantitative platform for verifying predicted splicing isoforms [110] [111]. This protocol details the application of this integrated technique to validate candidate AS events identified from bulk RNA-seq data, thereby bridging computational discovery with molecular confirmation.

Quantitative Performance Comparison of Detection Methods

The selection of a validation platform is guided by key performance metrics including sensitivity, specificity, reproducibility, and multiplexing capacity. The following table summarizes quantitative data comparing CE-based methods to other common nucleic acid analysis techniques, underscoring the advantages of the RT-PCR-CE workflow.

Table 1: Quantitative Comparison of Nucleic Acid Detection and Analysis Methods

Method	Target/Application	Key Performance Metrics	Advantages for Splicing Validation	Primary Reference
Capillary Electrophoresis with Laser-Induced Fluorescence (CE-LIF)	Quantification of circulating DNA in serum.	Intra-assay variability: 4.19%; Inter-assay variability: 6.91% [112] [113]. Lower variability than comparable real-time PCR.	High precision and reproducibility for fragment sizing and quantification.	Sang et al., 2006 [112] [113]
Real-Time Quantitative PCR (qPCR)	Quantification of circulating DNA in serum.	Intra-assay variability: 7.3%; Inter-assay variability: 14.92% [112] [113].	High sensitivity, but lower reproducibility for absolute quantification vs. CE.	Sang et al., 2006 [112] [113]
Multiplex PCR-Capillary Electrophoresis (MPCE)	HPV DNA genotyping in cervical cancer tissue.	Overall positivity rate: 97.9%; Superior detection of multiple infections vs. PCR-RDB (p=0.002) [111].	Excellent multiplexing capacity and specificity for distinguishing multiple targets/sizes in one reaction.	Qin et al., 2024 [111]
PCR-Reverse Dot Blot (PCR-RDB)	HPV DNA genotyping (targeting L1 gene).	Overall positivity rate: 92.9% [111]. Kappa agreement with MPCE: 0.390 (fair) [111].	Lower multiplexing capability and potential for false negatives with fragmented DNA.	Qin et al., 2024 [111]
Capillary Electrophoresis-based Multiplex RT-PCR	Detection of 13 upper respiratory tract pathogens.	Minimum detectable concentration: ≤ 0.5 copies/µL for 10/13 targets [110].	Exceptional sensitivity and broad, multiplexed detection suitable for analyzing multiple splice variants.	Gong et al., 2025 [110]

Integrated Protocol for Validating Alternative Splicing Events

This protocol provides a detailed workflow from computational identification of AS events to their experimental validation.

Computational Identification of Candidate Splicing Events from Bulk RNA-seq

• Input Data: Processed bulk RNA-seq data (BAM files) from control and experimental conditions, with at least four biological replicates per condition recommended for robust statistical power [62].
• Splicing Analysis:
  • Utilize a specialized software tool such as rMATS (version 3.2.5 or higher) to analyze aligned reads [62].
  • Specify the five major splicing event types: Skipped Exon (SE), Alternative 5' Splice Site (A5SS), Alternative 3' Splice Site (A3SS), Mutually Exclusive Exons (MXE), and Retained Intron (RI) [2].
  • For each event, the software calculates the Percent Spliced-In (PSI or Ψ) value, representing the proportion of transcripts including the alternative cassette.
  • Perform statistical testing (e.g., likelihood ratio test in rMATS) to identify differential alternative splicing (DAS) events between conditions. Filter results based on a False Discovery Rate (FDR) < 0.05 and an absolute ΔPSI > 0.10 (or a user-defined threshold relevant to the biological context) [62].
• Candidate Selection: Prioritize 3-10 high-confidence DAS events for validation, considering factors like PSI confidence intervals, gene functional relevance, and amplicon size suitability for subsequent PCR.

RNA Isolation and Quality Control

• Extraction: Isolate total RNA from the same biological samples used for RNA-seq using a column-based kit (e.g., RNeasy Plus). Include a DNase I digestion step to eliminate genomic DNA contamination [62].
• Quality Assessment:
  • Quantify RNA using a fluorometric method.
  • Assess integrity via an automated electrophoresis system (e.g., Agilent Bioanalyzer). Proceed only with samples possessing an RNA Integrity Number (RIN) > 7.0 [19].

Primer Design for Splice Isoform Discrimination

• Strategy: Design one pair of primers per splicing event. Flanking primers should bind in constitutive exons upstream and downstream of the alternative region. This strategy ensures amplification of all isoforms, yielding products of different sizes that can be resolved by CE [62].
• Design Criteria:
  • Amplicon Size: The difference in length between inclusion and exclusion products should be >15 base pairs for clear resolution by CE. Optimal total product size is 100-500 bp.
  • Specificity: Verify primer sequence uniqueness using BLAST against the appropriate genome.
  • Properties: Primer length: 18-22 bp; Tm: 58-62°C; avoid secondary structures. Include a universal fluorescent label (e.g., 6-FAM) on the 5’ end of the forward primer for CE detection.

Multiplex RT-PCR Amplification

• Reverse Transcription: Synthesize cDNA from 500 ng of total RNA using a reverse transcriptase with high processivity (e.g., SuperScript IV). Use a mixture of oligo(dT) and random hexamer primers to ensure full-length and unbiased representation.
• PCR Setup:
  • For each sample, set up a multiplex PCR reaction containing primer pairs for 2-4 different splicing events, provided their amplicon sizes are distinct and do not overlap.
  • Use a high-fidelity, hot-start DNA polymerase to minimize nonspecific amplification.
  • Include a positive control (cDNA with a known splicing pattern) and a no-template control (NTC).
• Thermocycling Parameters:
  • Initial Denaturation: 95°C for 2 min.
  • 35 Cycles: Denature at 95°C for 30 sec, Anneal at 60°C for 30 sec, Extend at 72°C for 45 sec.
  • Final Extension: 72°C for 5 min.

Capillary Electrophoresis Analysis and Quantification

• Sample Preparation: Dilute 1 µL of the multiplex PCR product in 9 µL of Hi-Di formamide containing a proprietary internal size standard (e.g., GeneScan 500 LIZ).
• Instrument Run: Perform electrophoresis on a genetic analyzer (e.g., Applied Biosystems 3500 Series). Use a performance-optimized polymer (e.g., POP-7) and the appropriate run module for fragment analysis.
• Data Interpretation:
  • Analyze raw data using fragment analysis software (e.g., GeneMapper).
  • The software will identify peaks corresponding to the fluorescently labeled amplicons, assigning a size (in base pairs) to each based on the internal standard.
  • For each splicing event, identify the peaks corresponding to the inclusion isoform and the exclusion isoform.
• Calculation of Experimental PSI (Ψ_exp):
  • For each event, calculate the peak height or peak area for the inclusion and exclusion products.
  • Calculate: Ψexp = (Signalinclusion) / (Signalinclusion + Signalexclusion)
  • This experimentally derived Ψexp is directly comparable to the computational PSI (Ψcomp) estimated from RNA-seq data by tools like rMATS.

Statistical Validation

• Correlation Analysis: Perform linear regression or calculate the Pearson correlation coefficient between the Ψcomp (from RNA-seq) and Ψexp (from RT-PCR-CE) values across all validated events and biological replicates. A strong positive correlation (e.g., R² > 0.85) confirms the accuracy of the computational pipeline.
• Differential Splicing Confirmation: Use a two-tailed Student's t-test to compare the Ψ_exp values between control and experimental groups for each event. Successful validation is achieved if the statistical significance (p < 0.05) and the direction of change (ΔΨ) match the computational prediction.

Workflow and Strategy Visualization

Workflow for Splicing Validation

Primer Strategy for Detecting Splice Isoforms

Research Reagent Solutions

Table 2: Essential Reagents and Materials for RT-PCR-CE Splicing Validation

Category	Reagent/Kit	Function in Protocol	Key Features
RNA Isolation	RNeasy Plus Micro/Mini Kit (Qiagen)	Isolation of high-quality total RNA from cells or tissue.	Includes gDNA eliminator column for effective genomic DNA removal [62].
cDNA Synthesis	SuperScript IV Reverse Transcriptase (Thermo Fisher)	Generation of first-strand cDNA from RNA templates.	High thermal stability and processivity for robust cDNA yield, especially for long or structured transcripts.
PCR Amplification	KAPA HiFi HotStart ReadyMix (Roche)	High-fidelity amplification of specific splice isoforms.	Hot-start enzyme minimizes nonspecific priming; high fidelity ensures accurate amplification.
Fluorescent Labeling	6-FAM Amidite (Integrated DNA Technologies)	Chemical synthesis of 5’-labeled oligonucleotide primers.	Fluorescent dye (6-Carboxyfluorescein) is a standard, stable label compatible with CE laser detectors.
Capillary Electrophoresis	GeneScan 500 LIZ Dye Size Standard (Thermo Fisher)	Internal lane standard for precise fragment sizing during CE.	Contains fragments of known length labeled with a distinct dye (LIZ) for accurate base pair calling.
	Hi-Di Formamide (Thermo Fisher)	Denaturing agent for sample preparation prior to CE.	Denatures PCR products into single strands for accurate size-based separation.
Analysis Software	GeneMapper Software v6 (Thermo Fisher)	Fragment analysis and peak calling from CE raw data.	Automatically assigns sizes, identifies peaks, and calculates peak height/area for quantitative analysis.
Multiplex Assay Design	SureX HPV Genotyping Kit (Ningbo HEALTH Gene Tech) [Example]	Commercial example of a multiplex PCR-CE assay for pathogen genotyping [111].	Demonstrates the principle of designing a single-tube PCR to amplify multiple targets with distinct sizes for CE resolution.

Digital Droplet PCR (ddPCR) for Absolute Quantification of Splice Isoforms

Abstract Within the framework of a thesis on alternative splicing (AS) analysis via bulk RNA sequencing (RNA-seq), this article details the application of droplet digital PCR (ddPCR) for the absolute quantification of mRNA splice isoforms. RNA-seq enables genome-wide discovery of splicing variations but often requires precise, targeted validation of specific isoforms. ddPCR addresses this need by providing absolute copy number quantification without standard curves, offering high precision and sensitivity for low-abundance targets. These Application Notes and Protocols synthesize current methodologies, from RNA-seq-informed primer design using tools like IsoPrimer to optimized ddPCR workflows for splice junctions. We present detailed experimental protocols, highlight applications in basic research and therapeutic development (e.g., exon-skipping therapies for Duchenne Muscular Dystrophy), and provide a framework for integrating ddPCR as an essential validation tool within a comprehensive AS research pipeline [114] [115] [116].

Alternative splicing is a fundamental mechanism for enhancing proteomic diversity, with the majority of human multi-exon genes producing multiple mRNA isoforms [114] [117]. Disruptions in splicing are implicated in numerous diseases, including cancer, neurodegenerative disorders, and muscular dystrophies [15] [116]. Bulk RNA-seq has become the cornerstone for genome-wide profiling of AS, enabling the discovery of splicing events, such as local splicing variations (LSVs) and changes in percent spliced in (Ψ) [15]. However, translating RNA-seq findings into biologically validated conclusions presents significant challenges. Short-read sequencing struggles to resolve full-length isoforms, and computational quantification involves inherent uncertainty in assigning reads to highly similar transcripts [117] [18].

Consequently, targeted validation of key splicing events is a critical step. While reverse transcription quantitative PCR (RT-qPCR) has been the traditional gold standard, droplet digital PCR (ddPCR) offers a paradigm shift for absolute quantification [114]. ddPCR partitions a sample into thousands of nanoliter-sized droplets, performing an endpoint PCR in each. The binary readout (positive/negative for amplification) allows for absolute quantification of target molecules via Poisson statistics, eliminating the need for external standard curves and reducing sensitivity to PCR efficiency variations [114] [116]. This is particularly advantageous for splice isoforms, where distinguishing between sequences that differ by only a few nucleotides requires exceptional specificity and where transcripts of interest may be expressed at very low levels [118] [119].

This document frames ddPCR methodology within the broader context of an AS research thesis that begins with bulk RNA-seq discovery. We provide detailed application notes and wet-lab protocols to enable researchers to seamlessly transition from bioinformatic identification of splicing events to their precise molecular validation.

Core Principle: ddPCR for Splice Junction Quantification

The absolute quantification of splice isoforms by ddPCR hinges on the precise targeting of unique exon-exon junctions. Unlike gene expression assays that target constitutive exons, isoform-specific assays must be designed to span the junction created by splicing, ensuring amplification of only the intended isoform [115] [116].

• Assay Design: The ideal target is a sequence spanning the junction where two exons connect. TaqMan probes with the fluorophore/quencher system are typically used, with the probe designed to hybridize across the junction itself. This ensures that fluorescence occurs only if the exact spliced product is present [116]. For highly similar isoforms, primer placement in unique constitutive or alternative exons is critical [114].
• Digital Partitioning and Absolute Quantification: In ddPCR, the reaction mixture is partitioned into approximately 20,000 droplets. Following PCR amplification, each droplet is analyzed for fluorescence. Droplets containing the target sequence (positive) are counted against those without it (negative) [116]. The absolute concentration of the target molecule in the original sample, measured in copies per microliter, is calculated using Poisson distribution statistics to account for droplets containing more than one template molecule. This direct counting method provides an absolute measure, in contrast to the relative quantification (Cq values) of qPCR [114] [118].

Table 1: Comparison of Methods for Splice Isoform Quantification.

Feature	Bulk RNA-seq (Short-Read)	RT-qPCR	Droplet Digital PCR (ddPCR)
Quantification Type	Relative (estimates Ψ or TPM)	Relative (requires standard curve/reference genes)	Absolute (copies/μL)
Throughput	High (genome-wide)	Low to medium (targeted)	Low to medium (targeted)
Precision for Low Abundance	Low (depends on coverage)	Moderate (Cq values become unreliable >35)	High (binary counting)
Standard Curve Required	No	Yes	No
Primary Role in AS Thesis	Discovery & hypothesis generation	Traditional validation	High-precision validation & absolute quantification
Key Advantage	Unbiased, global profiling	Well-established, cost-effective	High sensitivity, precision, no calibration needed
Key Limitation	Indirect inference of isoforms	Relative measure, efficiency-dependent	Higher cost per sample, lower throughput

Integrated Workflow: From RNA-seq to ddPCR Validation

A robust AS analysis pipeline seamlessly connects high-throughput discovery with targeted validation. The following workflow integrates bulk RNA-seq and ddPCR.

Workflow for Splicing Analysis from RNA-seq to Validation

• RNA-seq Discovery Phase: Process bulk RNA-seq data through a splicing analysis tool (e.g., MAJIQ, rMATS) to identify genes with significant differential splicing between conditions [15]. Outputs include Local Splicing Variations (LSVs) and estimates of percent spliced in (Ψ) or change in Ψ (dΨ).
• Target Selection & Assay Design: Select high-priority splicing events for validation. Use RNA-seq data to inform the design of specific primers and probes. Tools like IsoPrimer automate this process by using the transcript quantification from RNA-seq (e.g., Kallisto output) to prioritize primer pairs that are common to the expressed isoforms of interest, ensuring accurate representation of the transcript pool [115].
• ddPCR Validation Phase: Execute the wet-lab protocol (detailed in Section 4) to obtain absolute copy numbers for each isoform in all relevant samples.
• Integrated Data Analysis: Correlate the absolute copy numbers from ddPCR with the relative Ψ values from RNA-seq. This step validates the computational predictions, provides biologically grounded quantitative data, and can identify potential biases in the RNA-seq analysis pipeline [116].

Detailed Experimental Protocols

Protocol I: Absolute Quantification of Isoforms in Low-Input Samples (e.g., Stem Cells)

This protocol, adapted for splice isoforms, is based on a microfluidic digital PCR (mdPCR) method for muscle stem cells and is directly applicable to chip-based ddPCR systems [118].

• Step 1: Cell Lysis and Reverse Transcription

  • Isolate cells of interest (e.g., by FACS). Lyse a defined number of cells (e.g., 500) directly in a PCR tube.
  • Perform reverse transcription using a cells-to-cDNA kit to generate first-strand cDNA.

• Step 2: cDNA Preparation and Exonuclease Digestion

  • Perform a limited-cycle PCR (e.g., 18 cycles) with a high-fidelity polymerase to generate double-stranded DNA (dsDNA) from the cDNA template.
  • Treat the product with Exonuclease I to remove excess single-stranded primers. Heat-inactivate the enzyme.

• Step 3: ddPCR Reaction Setup

  • Prepare the ddPCR reaction mix containing the digested dsDNA, TaqMan Fast Advanced Master Mix, and isoform-specific primer/probe sets. For a splice junction assay, the probe must be designed across the junction [116].
  • Follow manufacturer instructions to generate droplets (e.g., using a droplet generator). Transfer the emulsion to a 96-well PCR plate.

• Step 4: PCR Amplification and Droplet Reading

  • Run endpoint PCR with a thermal profile optimized for the TaqMan assay.
  • Load the plate into a droplet reader. The instrument will count the number of fluorescence-positive and negative droplets for each channel.

• Step 5: Data Analysis

  • Use the manufacturer's software (e.g., QuantaSoft) to apply a fluorescence amplitude threshold and calculate the absolute concentration (copies/μL) for each target isoform in the original reaction, based on Poisson statistics.
  • Normalize copies per isoform to the number of cells analyzed to obtain copies per cell [118].

Protocol II: High-Dimensional Multiplex ddPCR for Splice Variants

This advanced protocol outlines steps for multiplexing several splice variant targets in one reaction, based on work quantifying Hepatitis B virus splice variants [119].

• Step 1: Assay Design for Long Amplicons and Multiplexing

  • Design assays to generate amplicons of varying lengths (e.g., 0.3-2 kbp) covering different splice junctions (e.g., sp1, sp2, sp3) and a reference un-spliced transcript.
  • Use a common primer pair for all targets within a gene where possible, with unique internal probes labeled with different fluorophores. This controls for primer efficiency differences [119].
  • Empirically test >300 primer/probe concentration combinations to balance signal intensity and minimize cross-talk.

• Step 2: RNA Processing and Cleanup

  • Extract total RNA. Perform reverse transcription with a high-capacity kit.
  • Implement a post-RT cleanup step (e.g., column purification) to remove enzymes, salts, and primers that can inhibit downstream multiplex PCR.

• Step 3: Optimized Multiplex ddPCR Cycling

  • Prepare the ddPCR supermix with the optimized multiplex primer/probe set and template cDNA.
  • Generate droplets.
  • Use an increased cycle number (e.g., up to 80 cycles) to ensure efficient amplification of all long amplicons, especially the less abundant spliced variants [119].
  • Optimize denaturation and annealing temperatures to favor amplification of the spliced minority targets.

• Step 4: Data Acquisition and Analysis

  • Read the plate. Software will separate droplets into clusters based on their fluorescence signature in 2D (FAM/HEX) or multi-dimensional plots.
  • Manually verify cluster separation and set gates to count droplets for each specific splice variant and the reference. Calculate the ratio or absolute amount of each splice variant.

Table 2: The Scientist's Toolkit: Essential Reagents & Materials

Item	Function/Description	Example/Consideration
TaqMan Probe Assays	Isoform-specific detection. Probe must span the exon-exon junction.	Custom-designed; critical for specificity [116].
ddPCR Supermix for Probes	Optimized reaction mix for droplet generation and PCR.	Bio-Rad ddPCR Supermix for Probes (no dUTP).
Droplet Generator & Reader	Instrumentation for partitioning samples and reading fluorescence.	Bio-Rad QX200 system or equivalent.
RNA Extraction Kit	High-quality, DNase-treated RNA is essential.	Column-based kits (e.g., RNeasy) [114].
Reverse Transcription Kit	cDNA synthesis with high efficiency.	Use random hexamers and/or oligo-dT primers [114].
Nuclease-Free Water	All aqueous solutions must be nuclease-free.	Prevents RNA/DNA degradation.
Microfluidic Chips (mdPCR)	Alternative to droplets for partitioning.	Fluidigm 48.48 Dynamic Array IFCs [118].
Exonuclease I	Digests excess single-stranded primers post-pre-amplification.	Critical for low-input protocols to reduce background [118].

Key Applications in Research and Drug Development

Validation of RNA-seq Splicing Predictions

The primary application in an AS research thesis is the orthogonal validation of computational findings. A study comparing ddPCR and RT-qPCR for quantifying BRCA1 isoforms found both methods accurate, but highlighted ddPCR's independence from standard curves [114]. ddPCR provides a definitive molecular count for specific junctions identified by RNA-seq tools like MAJIQ, confirming the existence and abundance of predicted splicing changes [15].

Monitoring Exon-Skipping Therapies

ddPCR is the preferred method for quantifying the efficacy of antisense oligonucleotide (AON) therapies that induce therapeutic exon skipping, as in Duchenne Muscular Dystrophy (DMD) [116].

• Method: TaqMan assays are designed to target the novel exon-exon junction created by skipping (e.g., exon 50-52 junction when exon 51 is skipped) and the original wild-type junction.
• Advantage: ddPCR absolutely quantifies both the skipped and wild-type transcripts, accurately calculating the skipping percentage without the amplification bias that plagues nested PCR and the efficiency correction needed in qPCR [116]. This precision is critical for dose-response studies in pre-clinical and clinical development.

Profiling Viral or Cancer-Associated Splice Variants

Research on pathogenic splice variants benefits from ddPCR's sensitivity.

• Hepatitis B (HBV): High-dimensional multiplex ddPCR can simultaneously quantify multiple spliced HBV RNA isoforms (sp1, sp2, sp3) alongside the full-length pgRNA in serum and liver tissue, revealing their potential roles in infection and persistence [119].
• Cancer: ddPCR excels in detecting copy number variations (CNVs) of key exons in heterogeneous tumor samples. It has been shown to reclassify ambiguous MLPA results for BRCA1/2 exon-level deletions in advanced prostate cancer, offering more sensitive detection in samples with low tumor purity [120].

Data Analysis and Interpretation

• Threshold Setting: Proper analysis requires defining a fluorescence amplitude threshold to distinguish positive from negative droplets. For single assays, this can be set using no-template controls (NTCs). For multiplex assays, thresholds are set in 2D plots to separate distinct droplet clusters [119].
• Calculating Absolute Quantification: The software calculates the concentration (copies/μL) using the fraction of positive droplets (p) and the droplet volume (v): Concentration = −ln(1−p) / v. Results for different isoforms from the same sample can be expressed as a ratio (e.g., % exon skip = [skipped] / ([skipped]+[wild-type])) or as absolute copies per cell when normalized to cell count [118] [116].
• Integration with RNA-seq Data: Compare the ddPCR-derived absolute copy numbers or isoform ratios with the RNA-seq-derived Ψ values. A strong positive correlation validates the RNA-seq analysis pipeline for that event. Discrepancies may indicate limitations in the short-read assembly of that particular isoform or biases in the computational model [117].

Logic of RNA-seq-Informed Primer Design with IsoPrimer

Within a thesis centered on alternative splicing discovery via bulk RNA-seq, droplet digital PCR (ddPCR) emerges as an indispensable tool for closing the loop between computational prediction and biological reality. Its capacity for absolute, precise, and sensitive quantification of specific splice junctions provides a level of validation that relative quantification methods cannot match. The protocols and applications detailed herein—from low-input stem cell analysis to multiplexed viral variant profiling—provide a comprehensive roadmap for incorporating ddPCR into the AS research workflow. As splicing-directed therapies advance, the role of ddPCR in preclinical and clinical pharmacodynamics will only grow, solidifying its position as a cornerstone technology for precise molecular measurement in both basic and translational research [118] [116] [119].

Comparing Short-Read and Long-Read Sequencing Platforms for Splicing Analysis

In the context of a broader thesis on alternative splicing analysis using bulk RNA sequencing, the choice of sequencing platform is a fundamental decision that directly impacts the discovery, quantification, and biological interpretation of transcript isoforms. Alternative splicing (AS) is a crucial regulatory mechanism, estimated to affect over 90% of multi-exon human genes, and its dysregulation is implicated in numerous diseases, including cancer and neurodegenerative disorders [49]. While short-read sequencing (e.g., Illumina) has been the workhorse for transcriptomics due to its high throughput and accuracy, it struggles to resolve full-length isoforms, often requiring complex computational deconvolution [121]. In contrast, long-read sequencing platforms (e.g., Oxford Nanopore Technologies, PacBio) generate reads that span complete mRNA transcripts, enabling the direct observation of splicing variants but historically at the cost of higher error rates and lower throughput [122]. This application note provides a detailed technical comparison, standardized protocols, and analysis frameworks to guide researchers in selecting and implementing the optimal platform for splicing-focused research, with an emphasis on generating robust, biologically actionable data for drug discovery and development.

Technical Platform Comparison for Splicing Analysis

The selection between short-read and long-read sequencing involves balancing multiple technical parameters that directly influence splicing outcome metrics.

Table: Comparative Performance of Sequencing Platforms for Splicing Analysis

Feature	Short-Read (Illumina)	Long-Read (Nanopore)	Long-Read (PacBio Iso-Seq)	Implication for Splicing Analysis
Typical Read Length	50-600 bp [122]	500 bp -> 2.3 Mb (commonly 10-30 kb) [122]	Up to 50 kb library inserts; average ~30 kb [122]	Long reads span full-length transcripts, eliminating assembly.
Raw Read Accuracy	High (>99.9%) [122]	Moderate (95-99%); improved with basecalling [8] [122]	High for CCS reads (>99% with sufficient passes) [122]	Lower accuracy can obscure splice site boundaries. CCS requires sufficient sequencing depth.
Primary Splicing Output	Splice junction counts, exon inclusion levels (PSI) [49] [15]	Direct full-length isoform sequences [123] [121]	Direct full-length isoform sequences (CCS reads) [57]	Long reads provide immediate isoform structure; short reads infer isoforms statistically.
Novel Isoform Discovery	Limited; relies on junction detection & assembly [15]	High; direct sequencing reveals unannotated isoforms [123] [57]	High; CCS approach yields high-quality novel isoforms [57]	Long-read platforms dramatically expand the catalog of splicing variants.
Quantitative Performance	Excellent for gene-level and junction-level expression [121]	Good; correlates with short-read but with more technical noise [121]	Good; robust for identifying major isoforms [121]	Short-reads remain the gold standard for precise quantification of expression.
Typical Throughput	Very High (Billions of reads/run)	High (Millions to billions of reads/run) [121]	Moderate (Millions of reads/run)	Throughput dictates cost and depth of isoform detection per sample.
RNA Input Requirement	Low (ng scale)	Varies by protocol (Direct RNA: µg; PCR-cDNA: ng scale) [121]	Moderate to High (µg scale for optimal library prep)	Influences feasibility for precious clinical or single-cell samples.
Detects RNA Modifications	No (requires specialized protocols)	Yes, natively in Direct RNA mode [121]	Limited (detects some via kinetic analysis)	Nanopore Direct RNA uniquely allows concurrent m6A profiling and isoform analysis [121].

Table: Key Metrics from Comparative Splicing Studies

Study Reference	Platforms Compared	Key Splicing-Specific Finding	Quantitative Outcome
Mouse Retina Profiling [123]	Illumina short-read vs. ONT long-read (single-cell)	Long-reads excel at precise transcript isoform identification and discover many novel isoforms.	Identified 44,325 transcript isoforms; 38% were previously uncharacterized, 17% cell-subclass specific [123].
SG-NEx Project (7 cell lines) [121]	Illumina, ONT (3 protocols), PacBio Iso-Seq	Long-read sequencing more robustly identifies major transcript isoforms compared to short-read.	Long-read protocols provided direct isoform resolution, reducing deconvolution ambiguity seen in short-read data [121].
Chicken Myogenesis [57]	PacBio Iso-Seq & Illumina short-read	Integrated analysis provides a comprehensive view, with Iso-Seq greatly expanding known isoform diversity.	Generated a master transcriptome of 61,146 full-length isoforms; over 52% were novel (NIC/NNC categories) [57].
MAJIQ v2 Benchmark [15]	Analysis of heterogeneous short-read datasets	Capturing unannotated (de novo) splice junctions increases detected differential splicing events by >30%.	Highlights the importance of computational methods that account for novel splicing in short-read data.

Detailed Experimental Protocols

Protocol 1: Standard Bulk Short-Read RNA-seq for Splicing Analysis

This protocol is optimized for the accurate quantification of splice junctions and differential splicing event (DSE) detection from Illumina platforms [49] [124].

• Step 1: RNA Isolation and Quality Control (QC)

  • Isolate total RNA using a guanidinium thiocyanate-phenol method or column-based kits. For formalin-fixed paraffin-embedded (FFPE) samples, use specialized extraction protocols.
  • Assess RNA integrity using an Agilent Bioanalyzer or TapeStation. A minimum RNA Integrity Number (RIN) of 8.0 is recommended for standard poly-A selection protocols. For ribo-depletion protocols, a RIN > 7.0 may be acceptable [124]. Low RIN leads to 3' bias and uneven coverage, compromising splice junction detection.

• Step 2: Library Preparation

  • Poly-A Selection: Use oligo(dT) beads to enrich for polyadenylated mRNA. This is standard for gene expression and splicing but excludes non-polyadenylated RNAs.
  • Ribo-Depletion: Use probe-based methods to remove ribosomal RNA (rRNA). This retains non-coding RNAs and pre-mRNA, providing information on intron retention and expanding splicing analysis [124].
  • Strand-Specific Protocol: Employ dUTP or adaptor-ligation methods to preserve strand information. This is critical for accurate annotation of antisense transcripts and overlapping genes.
  • Fragment RNA or cDNA to a target size of 200-500 bp. Perform reverse transcription, second-strand synthesis, end repair, A-tailing, and adapter ligation. Include unique dual indices (UDIs) for sample multiplexing.

• Step 3: Sequencing

  • Sequence on an Illumina NovaSeq, HiSeq, or NextSeq platform.
  • Use paired-end sequencing (2x 100 bp or 2x 150 bp). The paired-end information is vital for accurate alignment across splice junctions.
  • Recommended Depth: For standard differential gene expression, 30-40 million aligned reads per sample is sufficient. For splicing analysis, aim for 50-100 million aligned reads per sample to ensure sufficient coverage across splice junctions for robust PSI (Percent Spliced In) calculation [15].

Protocol 2: Bulk Long-Read RNA-seq for Isoform Resolution

This protocol outlines the process for full-length isoform sequencing using Oxford Nanopore Technologies (ONT), highlighting the choice between cDNA and direct RNA methods [121].

• Step 1: RNA QC and Input Requirement

  • Follow the same stringent QC as in Protocol 1. Input requirements vary:
    • PCR-cDNA Protocol: Requires the least input (often as low as 10-50 ng total RNA). Most accessible but introduces amplification bias.
    • Direct cDNA Protocol: Requires ~500 ng - 1 µg of high-quality RNA. Avoids PCR but includes reverse transcription.
    • Direct RNA Protocol: Requires the most input (≥1 µg of poly-A RNA). Sequences native RNA, preserving base modifications (e.g., m6A) [121].

• Step 2: Library Preparation (ONT PCR-cDNA)

  • Reverse Transcription: Use a poly(T) primer with a known 5' sequence tag to generate first-strand cDNA.
  • PCR Amplification: Amplify cDNA with a primer matching the 5' tag and a poly(T) primer. This step incorporates ONT-specific adapters.
  • Purification and Size Selection: Clean PCR product. Optional size selection (e.g., with SPRI beads) can remove very short fragments.
  • Adapter Ligation: Ligate ONT sequencing adapters (RAP or LSK) to the double-stranded cDNA.
  • Prime and Load: Add sequencing primer, load the library onto a FLO-PRO or similar flow cell.

• Step 3: Sequencing and Basecalling

  • Sequence on a MinION, GridION, or PromethION device.
  • Perform real-time basecalling using Guppy (ONT's production basecaller). Select the appropriate high-accuracy (HAC) model matching your flow cell and kit chemistry.
  • For retrospective, more accurate basecalling, use the super-accuracy (SUP) model, which is more computationally intensive [122].
  • Recommended Output: Target a minimum of 5-10 million aligned reads per sample for isoform discovery and quantification in bulk samples. Include spike-in controls (e.g., Sequin, SIRV, ERCC) in the library to assess quantitative accuracy and technical performance [121].

Protocol 3: Integrated Short- and Long-Read Analysis for Validation

This hybrid protocol leverages the strengths of both platforms for high-confidence isoform discovery and validation [123] [57].

• Parallel Sequencing: Prepare and sequence the same RNA sample(s) using both the standard short-read protocol (Protocol 1) and a long-read protocol (Protocol 2).
• Long-Read Guided Annotation:
  • Process long-reads through a pipeline (e.g., minimap2 alignment -> TALON or FLAIR for isoform clustering -> SQANTI3 for quality control and classification) [57].
  • Use the high-confidence, full-length isoforms discovered by long-reads to create a sample-specific transcriptome annotation (GTF file).
• Short-Read Quantification Against Custom Annotation:
  • Align short-reads to the reference genome using a splice-aware aligner (e.g., STAR, HISAT2).
  • Perform transcript-level quantification (e.g., using Salmon or kallisto in alignment-free mode) using the custom annotation generated from long-read data. This increases the accuracy of short-read-based isoform quantification.
• Validation and Consolidation:
  • Compare PSI values for splicing events calculated from short-read junction counts (e.g., using MAJIQ) with isoform abundances derived from long-reads.
  • Consolidate findings, using long-read data to confirm the structure of differentially spliced isoforms identified from short-read data.

Bioinformatics Analysis Pipelines

Table: Selected Computational Tools for Splicing Analysis

Tool Name	Best For Platform	Analysis Type	Key Principle	Considerations
DEXSeq [49]	Short-read	Differential Exon Usage (DEU)	Fits a negative binomial GLM to exon count data, testing for changes in exon usage independent of gene expression.	Robust and widely used; works on exon-level counts; requires good annotation.
rMATS [49]	Short-read	Differential Splicing Events (DSE)	Models splicing event types (SE, RI, etc.) using a hierarchical Bayesian model to calculate posterior inclusion difference.	Event-based; intuitive output; can handle replicates.
MAJIQ v2 [15]	Short-read (large/heterog. datasets)	Local Splicing Variations (LSVs)	Builds a splice graph, quantifies LSVs (ΔPSI) using Bayesian inference, and includes non-parametric tests (MAJIQ HET) for heterogeneous groups.	Powerful for complex and de novo splicing; scalable for 1000s of samples; includes VOILA visualization.
FLAIR	Long-read (Nanopore/PacBio)	Isoform Discovery & Quantification	Corrects, aligns, and collapses reads to define isoforms, then quantifies them in samples.	Designed for noisy long-reads; effective for identifying novel isoforms.
SQANTI3	Long-read (Nanopare/PacBio)	Isoform QC & Classification	Characterizes isoforms relative to reference annotation (FSM, ISM, NIC, NNC), filters artifacts, and provides structural classification.	Essential QC step before downstream analysis of long-read-derived isoforms.
Longcell [8]	Long-read (Single-cell)	Single-Cell Isoform Quantification	Recovers UMIs and cell barcodes from noisy reads, performs UMI-based error correction, and quantifies isoforms per cell.	Addresses key challenge of UMI "scattering" in Nanopore scRNA-seq; enables single-cell splicing analysis.

The Scientist's Toolkit: Essential Reagents and Materials

Table: Key Research Reagent Solutions for Splicing Analysis

Item	Function	Recommended Use/Product Example
RNA Integrity Number (RIN) Analyzer	Assesses RNA degradation. Critical for predicting library complexity and coverage uniformity.	Agilent Bioanalyzer or TapeStation [124]. Requirement: RIN > 8 for optimal splicing libraries.
Poly-A Selection Beads	Enriches for polyadenylated mRNA from total RNA, reducing ribosomal RNA background.	NEBNext Poly(A) mRNA Magnetic Isolation Module or Illumina Stranded mRNA Prep kits.
Ribo-depletion Kits	Removes ribosomal RNA without biasing towards poly-A tail, retaining non-coding and nascent RNA.	Illumina Ribo-Zero Plus or QIAseq FastSelect. Crucial for studying intron retention [124].
Strand-Specific Library Prep Kit	Preserves the original orientation of the transcript during cDNA synthesis.	Illumina Stranded Total RNA Prep or NEBNext Ultra II Directional RNA Library Prep. Mandatory for accurate annotation.
Spike-In Control RNA Mixes	Added to samples before library prep to monitor technical performance, quantify absolute expression, and assess sensitivity.	ERCC RNA Spike-In Mix (Ambion) or Sequins [121]. SIRVs (Lexogen) are specially designed for isoform-level controls [121].
Long-Read Library Prep Kits	Converts RNA into a format compatible with long-read sequencers. Choice depends on goal and input.	ONT PCR-cDNA Kit (low input). ONT Direct RNA Kit (for modifications). PacBio Iso-Seq Kit [121] [57].
High-Accuracy Basecaller	Converts raw electrical (ONT) or fluorescence (PacBio) signals into nucleotide sequences.	Guppy (HAC/SUP models) for ONT [122]. ccs for PacBio HiFi read generation [122] [57].
Splicing Analysis Software	Detects and quantifies changes in splicing patterns from sequence data.	DEXSeq, rMATS (short-read) [49]. MAJIQ (complex short-read analysis) [15]. FLAIR, SQANTI3 (long-read) [57].

The landscape of splicing analysis is powerfully served by both short-read and long-read sequencing platforms, each with distinct and complementary strengths. Short-read sequencing remains the optimal choice for high-precision, cost-effective quantification of splicing dynamics across large sample cohorts, supported by a mature ecosystem of statistical tools like DEXSeq and rMATS [49]. Conversely, long-read sequencing is transformative for de novo isoform discovery, resolving complex splicing patterns, and validating transcript structures, as evidenced by studies where over 50% of identified isoforms were novel [123] [57]. For a thesis centered on bulk RNA sequencing in alternative splicing research, the most robust strategy involves an integrated approach. This entails using long-read data to build a sample-informed transcriptome annotation, which then guides the more quantitative analysis of deep short-read datasets. This hybrid methodology maximizes confidence in both the discovery and quantification of splicing variants, providing a comprehensive foundation for understanding their role in development and disease—a critical step toward identifying novel therapeutic targets.

Abstract The accurate identification of differential alternative splicing (AS) events from bulk RNA sequencing (RNA-seq) data is fundamental to understanding gene regulation in development and disease. With a plethora of computational tools available, systematic benchmarking is essential to guide tool selection and advance methodological standards. This application note synthesizes findings from recent comprehensive benchmarking studies to evaluate the accuracy, reproducibility, and scalability of AS detection tools. We present standardized experimental protocols for performance validation, provide a decision framework for tool selection based on analytical priorities, and outline a suite of essential research reagents and computational resources. The discussion is framed within the broader thesis that robust, standardized benchmarking is a critical prerequisite for generating biologically meaningful and translatable insights into splicing dysregulation.

Alternative splicing is a key mechanism for enhancing transcriptomic and proteomic diversity, with precise regulation being critical for cellular homeostasis and implicated in numerous diseases [101]. The advent of bulk RNA-seq has made genome-wide detection of differential splicing (DS) events a mainstream analytical task, prompting the development of dozens of bioinformatic tools [125]. However, this abundance of choice presents a significant challenge: studies have reported immense discrepancies in the outputs of different tools when analyzing the same dataset [125]. Without standardized evaluation, tool selection can become arbitrary, compromising the reproducibility and biological validity of research findings.

This article addresses the core challenge of tool evaluation within the context of bulk RNA-seq research. We define the three pillars of effective benchmarking—accuracy (correct identification of true splicing changes), reproducibility (consistent results across replicates and platforms), and scalability (efficient handling of large datasets). The following sections detail methodologies for rigorous benchmarking, synthesize quantitative performance data from independent studies, and provide practical protocols for researchers.

Methodologies for Comprehensive Benchmarking

Effective benchmarking requires a structured, multi-faceted approach that moves beyond simple comparisons. The methodologies below outline best practices derived from recent large-scale evaluations [125] [126].

Benchmarking accuracy necessitates a ground truth. Two primary data strategies are employed:

• Simulated Data: Tools like ASimulatoR introduce controlled, known differential splicing events into RNA-seq reads [126]. This provides exact truth sets for calculating precision and recall but may not fully capture real-world sequencing complexities.
• Curated Real Data: Using well-characterized experimental datasets (e.g., from model organisms or cell lines with validated splicing perturbations) offers biological authenticity. However, the complete ground truth is often incomplete, making accuracy metrics more approximate.

A robust benchmarking study should utilize both. For example, a 2023 benchmark of 21 tools used simulated data to measure core accuracy and real cancer datasets to assess biological relevance [125].

Core Performance Metrics

Tools should be evaluated across a matrix of quantitative and qualitative metrics [127] [128]:

• Accuracy Metrics: Precision (the fraction of predicted DS events that are true) and Recall (the fraction of true DS events that are detected) are fundamental. The F1-score, their harmonic mean, provides a single balanced metric.
• Reproducibility Metrics: Measured through the consistency of results across technical replicates, different sequencing depths, or by running the same tool with different reasonable parameter settings.
• Scalability & Efficiency: Runtime and memory usage under different sample sizes (e.g., n=3 vs. n=100 samples) are critical for practical application. Update frequency of underlying genome annotations also impacts long-term usability [127].
• Functional Usefulness: The ability to discover novel splicing events beyond annotated splice graphs and to provide interpretable output formats (e.g., Percent Spliced In (PSI) values) for downstream biological analysis are key differentiators [125].

Results from Systematic Benchmarking Studies

Synthesizing results from independent benchmarks reveals clear performance trends and trade-offs among popular tools.

Performance of Event-Based Differential Splicing Tools

A landmark 2023 study benchmarked 21 algorithms on simulated data. The results showed significant performance variation, identifying a subset of tools that consistently outperformed others [125].

Table 1: Performance Summary of Select Event-Based DS Detection Tools (Based on Simulated Data Benchmarking) [125]

Tool	Primary Methodology	Key Strength	Reported Limitation	Typical Use Case
SUPPA2	PSI calculation from transcript quantification; fast.	Very high speed, good precision.	Relies on accurate isoform quantification.	Large-scale screening where speed is critical.
rMATS	Statistical modeling of junction counts.	High sensitivity, robust statistical framework.	Longer runtimes on large datasets.	Standardized DS analysis with detailed statistical output.
LeafCutter	Intron-centric clustering; does not require annotation.	Discovers novel splicing events, good performance.	Focuses on intron excision events only.	Exploratory analysis in less-annotated genomes or conditions.
DARTS	Integrates RNA-seq data with Bayesian deep learning.	High accuracy by integrating multiple data types.	Computational complexity, requires more input.	When auxiliary data (like CLIP-seq) is available for training.
MAJIQ	Quantification and visualization of complex splicing variations.	Handles complex local splicing variations (LSVs) well.	Steeper learning curve, resource-intensive.	Detailed analysis of specific loci with complex splicing patterns.

Key Finding: No single tool universally outperformed all others across every metric. The study concluded that SUPPA, DARTS, rMATS, and LeafCutter were top performers in the event-based category [125]. Furthermore, combining results from two complementary high-performing tools (a consensus approach) was shown to improve overall reliability [125] [126].

The Impact of Analytical Strategy on Performance

Benchmarking reveals that performance is not solely a function of the tool but also of the analytical pipeline.

• Splice-Aware Mapping: The initial read alignment step is crucial. The DICAST benchmark found that aligners like STAR and HISAT2 offered the best balance of accuracy and runtime for downstream splicing analysis [126].
• Low-Expression Filtering: Filtering out very lowly expressed transcripts or junctions before DS analysis was demonstrated as a simple yet effective method to reduce false positive calls [125].
• Tool-Pair Combination: Aggregating results from two tools with different algorithmic foundations (e.g., one high-precision and one high-recall tool) can yield a more robust final set of DS events [125].

Application Notes and Protocols

Based on benchmarking evidence, we propose the following protocols to ensure accurate, reproducible, and scalable splicing analysis.

A Standardized Workflow for Differential Splicing Analysis

The diagram below outlines a consensus workflow integrating best practices for accuracy and reproducibility.

Diagram Title: Standardized Bulk RNA-seq Splicing Analysis Workflow

Protocol for Tool Selection and Consensus Building

Selecting the optimal tool requires aligning tool capabilities with research goals. The following decision framework aids in this process [125] [126].

Diagram Title: Decision Protocol for Splicing Tool Selection

The Scientist's Toolkit: Essential Research Reagent Solutions

A robust splicing analysis requires both biological and computational reagents. The following table details essential components.

Table 2: Research Reagent Solutions for Splicing Analysis

Item	Function / Purpose	Example / Note
High-Quality RNA Samples	Starting material for library prep. Integrity (RIN > 8) is critical for accurate splice junction detection.	Isolated from tissues/cells of interest using guanidinium thiocyanate-phenol-chloroform extraction.
Strand-Specific RNA-seq Library Prep Kit	Prepares cDNA libraries that preserve the original strand orientation of transcripts, crucial for accurate splice graph reconstruction.	Kits from Illumina, NEB, or Takara Bio.
Splice-Aware Aligner Software	Maps sequencing reads to the reference genome, correctly spanning introns and identifying novel junctions.	STAR [126], HISAT2 [126]. These are benchmarked for performance balance.
Reference Transcriptome	A set of annotated transcripts (splice variants) for a genome. Serves as the basis for quantification and event annotation.	ENSEMBL or GENCODE annotations. Must be matched to the genome build used for alignment.
Differential Splicing Tool Suite	Core software for statistical detection of splicing changes between conditions.	rMATS (sensitivity) [125], SUPPA2 (speed) [125], LeafCutter (novelty) [125].
Benchmarking Dataset (Simulated)	Provides ground truth for validating tool performance on your own system or for method development.	Generated using tools like ASimulatoR [126].
Consensus Analysis Scripts	Custom code to merge and filter results from multiple DS tools to generate a high-confidence call set.	Implement in Python or R, based on published strategies [125].

Systematic benchmarking has clarified the landscape of differential splicing analysis, demonstrating that careful tool selection and pipeline design are non-negotiable for accurate and reproducible research. The consensus among studies is that employing a multi-tool consensus strategy, anchored by high-performing tools like rMATS, SUPPA2, or LeafCutter and preceded by rigorous alignment and filtering, significantly enhances result reliability [125] [126].

Future developments must address key challenges: 1) Standardization: The field requires unified reporting formats and benchmark datasets to facilitate direct tool comparison [126]. 2) Scalability: As study sizes grow, tools must balance algorithmic sophistication with computational efficiency. 3) Functional Integration: The next generation of benchmarks should evaluate not just statistical detection but also the ability to predict the functional protein-level consequences of splicing changes. By adhering to rigorous benchmarking principles, researchers can ensure their findings on alternative splicing are a solid foundation for downstream biological discovery and therapeutic development.

Integrating Splicing Data with Genomic and Proteomic Findings

This document provides detailed application notes and protocols for the integration of splicing data with genomic and proteomic findings. The content is framed within the context of a broader thesis on alternative splicing analysis in bulk RNA sequencing research, with the aim of identifying candidate genes and biomarkers for complex traits and diseases. Research demonstrates that genetic variants affecting splicing can generate proteomic diversity and impact phenotypic traits, often independently of changes in overall gene expression [129]. This makes the integration of multi-omics data a critical step for moving from correlative observations to causal mechanistic insights.

The core challenge addressed here is the reliable identification of biologically and clinically relevant splicing variants. While high-throughput RNA sequencing (RNA-seq) has made genome-wide splicing analysis possible, distinguishing functional, protein-altering splicing events from transcriptional noise requires orthogonal validation [129] [130]. This protocol details a proteogenomic approach, where RNA-seq-derived splicing data is integrated with genomic variant information and mass spectrometry-based proteomics to confirm the translation and differential expression of specific protein isoforms [131].

Foundational Concepts and Quantitative Data Integration

The Role of Splicing Quantitative Trait Loci (sQTL/irQTL)

A key genomic component in splicing analysis is the identification of genetic variants that regulate alternative splicing. These are known as splicing quantitative trait loci (sQTL) or isoform ratio QTLs (irQTLs) [130]. Mapping these loci involves associating genetic variants with changes in the relative abundance of transcript isoforms.

Recent large-scale studies, such as those in the Framingham Heart Study, have quantified the extensive genetic architecture of isoform variation. The following table summarizes key quantitative findings from an irQTL mapping study in human whole blood [130].

Table 1: Summary of irQTL Mapping Findings from a Whole-Blood Study (n=2,622 for discovery) [130]

QTL Type	Definition	Number Identified	Key Features and Examples
Common cis-irQTL	Variant within ±1 Mb of the isoform’s transcription start site (MAF ≥ 0.01).	>1.1 million variants for 10,883 isoforms (4,971 genes).	20% had no significant association with overall gene expression, indicating isoform-specific regulation. Most significant sentinel: rs4841 for RPS14 (explains 96% of isoform ratio variance).
Common trans-irQTL	Variant >1 Mb from TSS or on different chromosome (P < 1.5 × 10⁻¹³).	1,870 sentinel variants for 1,084 isoforms (590 genes).	31% of pairs involved genes on chromosome 6, largely in the HLA region. Most significant pair: rs1458255 (chr4) regulating RPL9P9 (chr15).
Rare cis-irQTL	Variant within ±1 Mb of TSS (0.003 < MAF < 0.01).	2,327 variants for 2,467 isoforms (1,428 genes).	Some rare variants explained a larger fraction of variance than common variants for the same gene (e.g., rs1128271 explained 38% of variance for a POLR2G isoform).

Selecting Tools for Differential Splicing Analysis from RNA-Seq

The first analytical step is detecting differential splicing events from bulk RNA-seq data. Multiple bioinformatic tools exist, differing in their statistical approaches and the unit of analysis (e.g., exon, transcript, or splicing event). The following table categorizes a selection of widely used and maintained tools, based on recent benchmark reviews [49].

Table 2: Comparison of Selected Differential Splicing Analysis Tools for Short-Read RNA-Seq [49]

Tool	Statistical Family	Primary Level of Analysis	Strengths and Notes
DEXSeq	Parametric (Negative Binomial GLM)	Exon	Tests for differential exon usage independent of gene expression. Robust and widely used; based on DESeq2 framework.
rMATS	Parametric (Bayesian)	Splicing Event	Quantifies five main event types (SE, RI, MXE, A5SS, A3SS). Handles replicates and provides false discovery rates.
SUPPA2	Non-parametric	Transcript/Event	Fast, alignment-free method using transcript-level quantification. Useful for large-scale screening.
LeafCutter	Non-parametric	Intron Cluster	Identifies differential intron splicing from junction reads without needing isoform annotation.
MAJIQ	Probabilistic	Splicing Junction/Variant	Models splice junction uncertainty to detect local splicing variations with confidence intervals.

Recommendation: For researchers new to the field, DEXSeq (for exon-level analysis) and rMATS (for event-level analysis) are recommended due to their high citation frequency, continued developer support, and comprehensive documentation [49].

Detailed Experimental Protocol: A Proteogenomic Integration Pipeline

The following protocol, adapted from the Splicify pipeline, provides a step-by-step methodology for integrating RNA-seq, genomic, and mass spectrometry data to detect differentially expressed protein isoforms resulting from alternative splicing [131].

Protocol Title: Detection of Differentially Expressed Protein Isoforms via RNA-Seq and Mass Spectrometry Integration

Objective: To identify and confirm disease- or condition-specific protein isoforms originating from alternative splicing by integrating bulk RNA-seq data with tandem mass spectrometry (MS/MS) data.

Part A: RNA-Seq Analysis and Isoform-Specific Database Generation

• RNA Sequencing & Quality Control:

  • Perform standard bulk RNA-seq library preparation and sequencing on your sample pairs (e.g., disease vs. control, treated vs. untreated).
  • Assess raw read quality using tools like FastQC [84] and perform adapter trimming/quality filtering with Trimmomatic [131] or Trim Galore [84].

• Read Alignment and Splicing Analysis:

  • Align cleaned reads to the reference genome using a splice-aware aligner such as STAR [131].
  • Perform differential splicing analysis using a tool like rMATS [131] to generate a list of statistically significant splicing events (e.g., skipped exons, retained introns) between conditions.

• Generate an Isoform-Specific Protein Sequence Database:

  • Extract the nucleotide sequences of the differentially spliced transcript variants identified by rMATS.
  • Perform in silico 3-frame translation of these variant sequences to predict all possible novel peptide sequences introduced or altered by the splicing event.
  • Combine these predicted novel peptide sequences with a standard human protein reference database (e.g., from UniProt [131]) to create a customized sample-specific protein FASTA database for mass spectrometry searching.

Part B: Mass Spectrometry Data Acquisition and Search

• Sample Preparation and Mass Spectrometry:

  • Prepare protein extracts from the same biological samples used for RNA-seq.
  • Perform tryptic digestion, peptide separation (typically via liquid chromatography), and analyze using tandem mass spectrometry (MS/MS) to generate fragmentation spectra.

• Database Search for Isoform-Specific Peptides:

  • Search the MS/MS spectra against the customized protein database generated in Step 3 using a search engine such as MaxQuant [131].
  • This search will identify peptides that map uniquely to the novel junctions or sequences created by alternative splicing. Three key peptide types are sought:
    • Junction-spanning peptides: Map across a novel exon-exon junction.
    • Inclusion-specific peptides: Map within a sequence exclusively present in an alternatively included exon or retained intron.
    • Frame-shift peptides: Result from a frameshift caused by the splicing event.

Part C: Integrative Data Analysis

• Validation and Quantification:

  • Filter the mass spectrometry results to retain only peptides that are unique to the splicing variant and have high-confidence identifications.
  • Quantify the abundance of these isoform-specific peptides across conditions using the intensity values from the MS data (e.g., MaxQuant output).
  • Correlate the peptide-level quantitative changes with the RNA-level differential splicing statistics from rMATS. A confirmed event shows significant differential splicing at the RNA level and corresponding differential detection/abundance of the isoform-specific peptide at the protein level.

• Prioritization and Functional Analysis:

  • Prioritize protein isoforms confirmed by multiple unique peptides or by peptides with high spectral counts.
  • Integrate genomic data (e.g., from whole genome sequencing) by checking if the splicing event is linked to a nearby genetic variant (e.g., an irQTL from Table 1 [130]). This can suggest a causal genetic mechanism.
  • Perform pathway analysis on the genes with confirmed differentially expressed protein isoforms to infer biological functions and potential therapeutic relevance.

Diagram: Splicify Proteogenomic Analysis Workflow [131]

Part D: Experimental Considerations and Troubleshooting

• Sample Matching: The most critical requirement is that the RNA and protein are extracted from the same biological specimen or from aliquots of an identical cell population to ensure valid correlation.
• Sensitivity: Mass spectrometry may not detect low-abundance protein isoforms. Deep proteomic coverage and sufficient sample loading are necessary.
• Peptide Specificity: Carefully validate that peptides claimed to be isoform-specific are not present in other proteins in the reference database. The customized database mitigates this risk.
• Replication: Biological replicates are essential at both the RNA-seq and proteomics stages to ensure statistical rigor.

The Scientist's Toolkit: Research Reagent Solutions

Table 3: Essential Materials and Tools for Splicing Data Integration Experiments

Item	Function/Description	Example/Supplier
Total RNA Isolation Kit	To extract high-quality, intact RNA from cells or tissues for RNA-seq. Integrity (RIN > 7-8) is crucial for accurate splicing analysis [19].	PicoPure RNA Isolation Kit (Thermo Fisher) [19].
Poly(A) mRNA Magnetic Beads	For enrichment of messenger RNA from total RNA prior to library preparation, reducing ribosomal RNA contamination.	NEBNext Poly(A) mRNA Magnetic Isolation Module [19].
Stranded RNA-seq Library Prep Kit	Prepares sequencing libraries from mRNA, preserving strand information which improves transcriptome annotation.	NEBNext Ultra II Directional RNA Library Prep Kit [19].
Splicing-Aware Aligner Software	Aligns RNA-seq reads to the genome, correctly mapping reads that span splice junctions. Essential for splicing detection.	STAR [131], HISAT2 [49].
Differential Splicing Analysis Software	Statistically tests for significant differences in splicing patterns between sample groups.	rMATS [131], DEXSeq [49].
Mass Spectrometry-Grade Trypsin	Proteolytic enzyme used to digest proteins into peptides for LC-MS/MS analysis.	Trypsin, sequencing grade (Promega).
Liquid Chromatography-Tandem Mass Spectrometer (LC-MS/MS)	Instrument platform for separating, ionizing, and fragmenting peptides to generate spectra for protein identification.	Orbitrap-series, Q1-TOF platforms.
Proteomics Search Engine Software	Matches experimental MS/MS spectra against a protein sequence database to identify peptides and proteins.	MaxQuant [131], Proteome Discoverer.
Proteogenomic Pipeline Scripts	Customizable computational workflow to integrate RNA-seq splicing calls with MS/MS identification.	Splicify (GitHub) [131].

Visualization of Integrative Analytical Pathways

The final confirmation of a splicing event's functional impact requires linking genetic regulation to transcript and protein output. The following diagram synthesizes the key steps in this integrative pathway, from genetic variant to phenotype [129] [130] [131].

Diagram: Integrative Path from Genetic Variant to Splicing-Mediated Phenotype

Alternative splicing (AS) is a fundamental regulatory mechanism, affecting over 90% of multi-exon human genes and serving as a critical source of proteomic diversity [8]. Its dysregulation is not a mere cellular artifact but a direct contributor to disease pathogenesis. An estimated 15–30% of all disease-causing mutations disrupt normal splicing patterns, contributing to a vast array of conditions including cancer, neurodegenerative disorders, and rare genetic diseases [13].

The convergence of high-throughput RNA sequencing (RNA-seq) and sophisticated computational biology has moved splicing analysis from basic research to the precipice of clinical utility. The thesis of this application note is that bulk RNA-seq data, when analyzed with splicing-aware pipelines, provides a powerful, reproducible, and actionable dataset for clinical decision-making. This is evidenced by the success of RNA-targeted therapies like nusinersen for spinal muscular atrophy and eteplirsen for Duchenne muscular dystrophy, which function by correcting disease-causing splicing errors [13]. This document provides detailed protocols and frameworks to systematically identify, quantify, and interpret splicing variations from bulk RNA-seq data, enabling researchers and drug developers to bridge the gap between transcriptomic discovery and clinical application.

Foundational Protocols: Bulk RNA-seq for Splicing Detection

Accurate splicing analysis begins with rigorous experimental design and data generation. Bulk RNA-seq remains the most robust and cost-effective method for profiling splicing patterns across patient cohorts or tissue samples [49].

Experimental Design & Wet-Lab Protocol

• Sample Collection & RNA Integrity: Use matched case-control or longitudinal sample sets. Preserve samples in RNA-stabilizing reagents immediately upon collection. Assess RNA quality using an Agilent Bioanalyzer or TapeStation; an RNA Integrity Number (RIN) > 8.0 is recommended for reliable splicing analysis [19].
• Library Preparation: Employ poly(A) selection to enrich for mature mRNA. Use strand-specific library preparation kits to accurately determine the origin of transcripts. While standard Illumina short-read sequencing (e.g., 150bp paired-end) is sufficient for most splicing event detection, consider long-read technologies (PacBio HiFi or Oxford Nanopore) for de novo isoform discovery in genes with complex architecture [2] [1].
• Sequencing Depth: For splicing-focused studies, prioritize depth over the number of samples. A minimum of 50-100 million reads per sample is recommended to ensure sufficient coverage across splice junctions, especially for lowly expressed genes [132].
• Controls & Replicates: Include positive control samples with known splicing variants if available. A minimum of three biological replicates per condition is essential for robust statistical power in differential splicing analysis [133].

Core Computational Pre-processing Workflow

The following workflow is optimized for accuracy and reproducibility. The logical flow from raw data to a quantified splicing event table is depicted in Figure 1.

Detailed Protocol Steps:

• Quality Control: Use FastQC to assess per-base sequence quality, adapter contamination, and GC content. Trim low-quality bases and adapter sequences using Trimmomatic or cutadapt [132].
• Splice-Aware Alignment: Align reads to the reference genome using a splice-aware aligner. STAR is highly recommended due to its speed and accuracy in detecting canonical and non-canonical splice junctions [49]. Use the --twopassMode Basic option for improved novel junction discovery.
• Quantification of Splicing Events: Quantify predefined splicing events (exon skipping, intron retention, alternative 5'/3' splice sites) using specialized tools. This step generates a count matrix for each splicing event (e.g., included-exon reads vs. skipped-exon reads) and a key metric: the Percent Spliced In (PSI or Ψ), which represents the proportion of transcripts containing a particular exon or splice site.

Computational Analysis: From Aligned Reads to Differential Splicing

The core challenge is to statistically distinguish biological splicing changes from technical noise. Multiple algorithmic approaches exist, as summarized in Table 1.

Table 1: Comparison of Major Differential Splicing Analysis Tools [49] [2]

Tool Name	Statistical Family	Primary Analysis Level	Key Strength	Consideration for Clinical Use
rMATS	Likelihood-based	Splicing Event (Exon-centric)	Models replicate variance; excellent for case-control studies.	Mature, widely used. Outputs PSI and P-values directly.
DEXSeq	Parametric (NB GLM)	Exon Bin	Detects differential exon usage independent of predefined event types.	High sensitivity but requires careful control of gene-level expression confounders.
LeafCutter	Non-parametric	Intron Clustering	Discovers novel splicing events from intron excision clusters without prior annotation.	Powerful for finding unannotated events in poorly annotated genomes or disease states.
SUPPA2	Fast-algorithmic	Transcript/Event	Uses transcript-level quantification (e.g., from Salmon) to rapidly compute PSI changes.	Extremely fast for large cohorts; depends on accuracy of transcript quantification.
JunctionSeq	Parametric (NB GLM)	Splice Junction	Focuses directly on read counts spanning splice junctions.	Complements exon-based methods; useful for alternative splice site detection.

Protocol: Differential Splicing with rMATS:

• Input: STAR-generated BAM files for all samples, along with a text file defining sample groups (e.g., Control vs. Disease).
• Run Command: Execute rMATS with parameters for paired-end reads and the desired statistical cutoff (e.g., FDR < 0.05).

• Output Interpretation: The primary output file (JC.raw.input.psi.txt) lists all tested splicing events, their PSI in each group, the difference in PSI (ΔPSI), and the associated False Discovery Rate (FDR). Filter for events with |ΔPSI| > 0.10 (10% change) and FDR < 0.05 as high-confidence differential splicing candidates [49].

Clinical Translation: Interpreting and Validating Splicing Variants

Identifying a differential splicing event is only the first step. Determining its clinical relevance requires a structured interpretation and validation pipeline (Figure 2).

Interpretation Protocol:

• Prioritization: Filter events based on effect size (ΔPSI), statistical significance (FDR), and gene function. Prioritize events in disease-associated genes and those predicted to cause nonsense-mediated decay (NMD) or major protein domain truncation [13].
• In Silico Pathogenicity Prediction: Use tools like SpliceAI or SPiP to assess whether the splicing alteration (or a genomic variant causing it) is predicted to be damaging. Integrate this evidence with the ACMG/ClinGen variant interpretation guidelines for splicing variants (PVS1 criterion) [13].
• Experimental Validation (Gold Standard):
  • RT-PCR & Capillary Electrophoresis: Design primers flanking the alternative exon. Amplify cDNA from original patient samples and controls. Resolve and quantify isoforms using capillary electrophoresis (e.g., Agilent Bioanalyzer) for high-resolution size discrimination and accurate quantification [1].
  • Minigene Splicing Assay: Clone the genomic region of interest (including exons and introns) into a splicing reporter vector. Transfert cells with the wild-type and mutant constructs, isolate RNA, and use RT-PCR to assess if the variant directly causes the aberrant splicing seen in the RNA-seq data [1].

Table 2: Key Research Reagent Solutions for Splicing Analysis

Item	Function	Example/Provider
RNA Stabilization Reagent	Preserves transcriptome integrity immediately upon sample collection, critical for accurate splicing profiles.	RNAlater (Thermo Fisher), PAXgene (PreAnalytiX)
Stranded mRNA Library Prep Kit	Generates sequencing libraries that preserve strand-of-origin information, crucial for accurate transcript annotation.	NEBNext Ultra II Directional (NEB), TruSeq Stranded mRNA (Illumina)
Splice-Aware Aligner Software	Aligns RNA-seq reads to a reference genome, accurately mapping reads across splice junctions.	STAR [132], HISAT2 [49]
Splicing Quantification Tool	Quantifies reads supporting different splicing events from aligned data (BAM files).	rMATS [49], DEXSeq [49]
In Silico Prediction Tool	Predicts the impact of genomic sequence changes on splicing.	SpliceAI [13], SPiP [13]
High-Fidelity Polymerase	Used for validation PCRs to minimize amplification errors and ensure accurate isoform detection.	Q5 (NEB), Phusion (Thermo Fisher) [1]
Capillary Electrophoresis System	Provides high-resolution, quantitative analysis of PCR products for splice isoform validation.	Agilent Bioanalyzer 2100, Fragment Analyzer [1]

Conclusion

Alternative splicing analysis with bulk RNA-Seq has matured into a powerful discipline essential for unraveling transcriptomic complexity in health and disease. A successful analysis hinges on a solid grasp of biological foundations, a carefully selected and executed methodological pipeline, rigorous optimization to mitigate false discoveries, and thorough validation of key findings. Future directions will be shaped by the integration of long-read sequencing to fully resolve isoform structures, the application of splicing biomarkers in clinical diagnostics, and the development of therapies that directly target splicing dysregulation. By adhering to the comprehensive framework outlined here, researchers can robustly capture the dynamic landscape of splicing, unlocking new avenues for biomarker discovery, therapeutic target identification, and the advancement of precision medicine.

References

• 1. Methods for Analyzing Alternative Splicing and Its Regulation ... [https://pmc.ncbi.nlm.nih.gov/articles/PMC12626760/]
• 2. Computational methods for alternative polyadenylation and ... [https://pmc.ncbi.nlm.nih.gov/articles/PMC12411624/]
• 3. The implications of alternative splicing regulation for ... [https://www.nature.com/articles/s41467-025-65339-1]
• 4. Predicting the structural impact of human alternative splicing [https://genomebiology.biomedcentral.com/articles/10.1186/s13059-025-03744-x]
• 5. Deciphering alternative splicing patterns during cell fate ... [https://pmc.ncbi.nlm.nih.gov/articles/PMC12150459/]
• 6. Long-read RNA sequencing: A transformative technology ... [https://www.sciencedirect.com/science/article/pii/S1525001624007524]
• 7. Computational methods for alternative polyadenylation and ... [https://www.rna-seqblog.com/computational-methods-for-alternative-polyadenylation-and-splicing-in-post-transcriptional-gene-regulation/]
• 8. Single cell and spatial alternative splicing analysis with ... [https://www.nature.com/articles/s41467-025-60902-2]
• 9. Alternative splicing [https://en.wikipedia.org/wiki/Alternative_splicing]
• 10. Alternative splicing in stem cells and development [https://pmc.ncbi.nlm.nih.gov/articles/PMC12137872/]
• 11. Mechanism of alternative splicing and its regulation - PMC [https://pmc.ncbi.nlm.nih.gov/articles/PMC4360811/]
• 12. RNA-Seq Analysis of an Antisense Sequence Optimized for ... [https://pmc.ncbi.nlm.nih.gov/articles/PMC5785776/]
• 13. Genome-wide functional annotation and interpretation of ... [https://www.nature.com/articles/s10038-025-01424-z]
• 14. Identification of alternative 5′/3′ splice sites based on ... [https://pmc.ncbi.nlm.nih.gov/articles/PMC1669764/]
• 15. RNA splicing analysis using heterogeneous and large ... [https://www.nature.com/articles/s41467-023-36585-y]
• 16. Incorporating exon–exon junction reads enhances differential ... [https://bmcbioinformatics.biomedcentral.com/articles/10.1186/s12859-025-06210-4]
• 17. Challenges in detecting and quantifying intron retention ... [https://www.sciencedirect.com/science/article/pii/S2001037019303721]
• 18. [Tutorial] Bulk RNA-seq DE analysis [https://informatics.fas.harvard.edu/resources/tutorials/differential-expression-analysis/]
• 19. A Beginner's Guide to Analysis of RNA Sequencing Data [https://pmc.ncbi.nlm.nih.gov/articles/PMC6096346/]
• 20. Alternative splicing as a biomarker and potential target for ... [https://pmc.ncbi.nlm.nih.gov/articles/PMC4648177/]
• 21. Identifying Splicing Regulatory Elements with de Bruijn Graphs [https://pmc.ncbi.nlm.nih.gov/articles/PMC4253301/]
• 22. A deep learning approach to identify gene targets of a ... [https://www.nature.com/articles/s41467-021-23663-2]
• 23. A General Definition and Nomenclature for Alternative ... [https://pmc.ncbi.nlm.nih.gov/articles/PMC2467475/]
• 24. Deciphering splicing heterogeneity at single-cell resolution ... [https://www.nature.com/articles/s41467-025-64517-5]
• 25. a resource for analysis and visualization of RNA-Seq data on ... [https://pmc.ncbi.nlm.nih.gov/articles/PMC3436850/]
• 26. Exploring RNA-Seq Data Analysis Through Visualization ... [https://www.mdpi.com/2306-5354/12/1/56]
• 27. Integrated analysis of single-cell RNA-seq and bulk RNA-seq ... [https://pmc.ncbi.nlm.nih.gov/articles/PMC12548165/]
• 28. Single-cell and bulk transcriptome profiling reveals RNA ... [https://www.nature.com/articles/s41598-025-14823-1]
• 29. A real-world multi-center RNA-seq benchmarking study ... [https://www.nature.com/articles/s41467-024-50420-y]
• 30. Integrative analysis of bulk and single-cell RNA ... [https://www.frontiersin.org/journals/medicine/articles/10.3389/fmed.2025.1552029/full]
• 31. RNA-Seq Experimental Design Guide for Drug Discovery [https://www.lexogen.com/blog/planning-for-success-a-strategic-design-guide-for-rna-seq-experiments-in-drug-discovery/]
• 32. Computational Comparison of Differential Splicing Tools ... [https://pmc.ncbi.nlm.nih.gov/articles/PMC11989494/]
• 33. Top Bioinformatics Tools for RNA Sequencing Data Analysis [https://riveraxe.com/top-bioinformatics-tools-for-rna-sequencing-data-analysis-best-value-for-your-research/]
• 34. Manual | HISAT2 [http://daehwankimlab.github.io/hisat2/manual/]
• 35. Mastering RNA-seq Read Alignment With STAR [https://www.reneshbedre.com/blog/star-aligner.html]
• 36. Evaluation of Seven Different RNA-Seq Alignment Tools Based on... [https://pmc.ncbi.nlm.nih.gov/articles/PMC7084517/]
• 37. EASTR: Identifying and eliminating systematic alignment ... [https://www.nature.com/articles/s41467-023-43017-4]
• 38. Optimal Dual RNA-Seq Mapping for Accurate Pathogen Detection in... [https://bio-protocol.org/en/bpdetail?id=5182&type=0]
• 39. Evaluation and comparison of computational tools for RNA-seq isoform quantification | BMC Genomics [https://link.springer.com/article/10.1186/s12864-017-4002-1]
• 40. Computational comparison of common event-based ... [https://bmcbioinformatics.biomedcentral.com/articles/10.1186/s12859-021-04263-9]
• 41. Comparative evaluation of full-length isoform quantification from RNA-Seq | BMC Bioinformatics [https://link.springer.com/article/10.1186/s12859-021-04198-1]
• 42. Systematic evaluation of differential splicing tools for RNA- ... [https://pubmed.ncbi.nlm.nih.gov/31802105/]
• 43. Evaluation of differential splicing tools for RNA-seq studies [https://www.rna-seqblog.com/evaluation-of-differential-splicing-tools-for-rna-seq-studies/]
• 44. Comprehensive assessment of mRNA isoform detection ... [https://www.nature.com/articles/s41467-024-48117-3]
• 45. Simulation-based benchmarking of isoform quantification in ... [https://genomebiology.biomedcentral.com/articles/10.1186/s13059-018-1571-5]
• 46. Isoform-level quantification for single-cell RNA sequencing [https://pmc.ncbi.nlm.nih.gov/articles/PMC8826380/]
• 47. Quantification of transcript isoforms at the single-cell level ... [https://www.nature.com/articles/s41467-025-61118-0]
• 48. Systematic evaluation of differential splicing tools for RNA-seq studies - PMC [https://pmc.ncbi.nlm.nih.gov/articles/PMC7711265/]
• 49. Selecting differential splicing methods - PubMed Central - NIH [https://pmc.ncbi.nlm.nih.gov/articles/PMC12308171/]
• 50. rMATS: Robust and flexible detection of differential ... [https://pmc.ncbi.nlm.nih.gov/articles/PMC4280593/]
• 51. LeafCutter: Annotation-free quantification of RNA splicing [https://davidaknowles.github.io/leafcutter/]
• 52. LeafCutterMD: an algorithm for outlier splicing detection in ... [https://pmc.ncbi.nlm.nih.gov/articles/PMC7750945/]
• 53. Tools To Do the Alternative Splicing Analysis [https://www.biostars.org/p/178032/]
• 54. Differential Alternative Splicing Analysis with rMATS [https://bioinformatics.ccr.cancer.gov/docs/btep-coding-club/CC2023/rmats/]
• 55. Comprehensive and scalable quantification of splicing differences with MntJULiP | Genome Biology | Full Text [https://genomebiology.biomedcentral.com/articles/10.1186/s13059-022-02767-y]
• 56. MAJIQ | Welcome [https://majiq.biociphers.org/]
• 57. Bulk and single‐cell alternative splicing analyses reveal roles ... [https://pmc.ncbi.nlm.nih.gov/articles/PMC10849790/]
• 58. PTBP1 depletion in mature astrocytes reveals distinct ... [https://elifesciences.org/articles/107683]
• 59. A spatial long-read approach at near-single-cell resolution ... [https://www.nature.com/articles/s41467-025-63301-9]
• 60. Discovery of RNA‐Targeting Small Molecules - PubMed Central [https://pmc.ncbi.nlm.nih.gov/articles/PMC12375692/]
• 61. Expression and splicing mediate distinct biological signals [https://pmc.ncbi.nlm.nih.gov/articles/PMC10588054/]
• 62. Bulk RNA-seq and alternative splicing analysis [https://bio-protocol.org/exchange/minidetail?id=8708131&type=30]
• 63. The hitchhikers' guide to RNA sequencing and functional ... [https://pmc.ncbi.nlm.nih.gov/articles/PMC9851315/]
• 64. Bulk RNA sequencing dataset of Claudin-low breast ... [https://www.nature.com/articles/s41597-025-06332-7]
• 65. Systematic review and meta-analysis of bulk RNAseq studies ... [https://pubmed.ncbi.nlm.nih.gov/40042520/]
• 66. A meta-analysis of bulk RNA-seq datasets identifies ... [https://www.nature.com/articles/s41598-024-75431-z]
• 67. Consistently processed RNA sequencing data from 50 ... [https://www.nature.com/articles/s41597-025-05376-z]
• 68. Protocol for unlocking alternative polyadenylation insights ... [https://www.sciencedirect.com/science/article/pii/S266616672500070X]
• 69. Current Progress and Future Perspectives of RNA-Based ... [https://pmc.ncbi.nlm.nih.gov/articles/PMC12153701/]
• 70. Power analysis and sample size estimation for RNA-Seq ... [https://pmc.ncbi.nlm.nih.gov/articles/PMC4201821/]
• 71. Sample Size Estimation for Detection of Splicing Events in ... [https://pmc.ncbi.nlm.nih.gov/articles/PMC5618549/]
• 72. Sample size, power and effect size revisited: simplified and ... [https://pmc.ncbi.nlm.nih.gov/articles/PMC7745163/]
• 73. How to do a power analysis to determine sample size [https://www.statsig.com/perspectives/power-analysis-determine-sample-size]
• 74. Sample size calculation based on exact test for assessing ... [https://bmcbioinformatics.biomedcentral.com/articles/10.1186/1471-2105-14-357]
• 75. Experimental design considerations | Introduction to RNA-seq ... [https://hbctraining.github.io/Intro-to-rnaseq-fasrc-salmon-flipped/lessons/02_experimental_planning_considerations.html]
• 76. Guidelines for bulk RNA Sequencing [https://genomecenter.tu-dresden.de/infohub/bulk_guidelines/bulkRNA-Seq.html]
• 77. Tumor Tissue Microarrays (TMAs): Importance, ... [https://blog.crownbio.com/tumor-tissue-microarrays-cancer-research]
• 78. 2025-2030 Global Automatic Tissue Microarray Analyzer ... [https://www.hnyresearch.com/report/2025-2030-Global-Automatic-Tissue-Microarray-Analyzer-Outlook-Market-Size-Share--Trends-Analysis-Report-By-Player-Type-Application-and-Region/1832528]
• 79. Dissecting Cellular Heterogeneity Using Single-Cell RNA ... [https://pmc.ncbi.nlm.nih.gov/articles/PMC6449718/]
• 80. Unveiling tissue heterogeneity through genomic interaction ... [https://www.sciencedirect.com/science/article/pii/S0002929725003556]
• 81. What is Bulk sequencing ( RNA - Bulk )? - Single Cell... RNA seq [https://www.scdiscoveries.com/support/what-is-bulk-rna-sequencing/]
• 82. Navigating single-cell RNA -sequencing: protocols , tools, databases... [https://pmc.ncbi.nlm.nih.gov/articles/PMC12085826/]
• 83. DiSC: a statistical tool for fast differential expression analysis ... [https://pubmed.ncbi.nlm.nih.gov/40444783/]
• 84. a versatile toolkit for analyzing and visualizing RNA-Seq data [https://bmcbioinformatics.biomedcentral.com/articles/10.1186/s12859-024-05761-2]
• 85. Multiple testing correction in RNA-Seq [https://www.biostars.org/p/9568684/]
• 86. What is multiple testing correction? [https://geneviatechnologies.com/blog/what-is-multiple-testing-correction/]
• 87. False Discovery Rate [https://www.publichealth.columbia.edu/research/population-health-methods/false-discovery-rate]
• 88. Benjamini-Hochberg correction: reducing false positives in ... [https://www.statsig.com/perspectives/benjamini-hochberg-false-positives]
• 89. How does multiple testing correction work? - PMC [https://pmc.ncbi.nlm.nih.gov/articles/PMC2907892/]
• 90. Beware of counter-intuitive levels of false discoveries in ... [https://pmc.ncbi.nlm.nih.gov/articles/PMC12359981/]
• 91. 转录组学的生物信息学视角：批量与单细胞RNA 测序分析的 ... [https://blog.sciencenet.cn/home.php?mod=space&uid=571917&do=blog&id=1512405]
• 92. RNA测序| RNA-Seq方法和工作流程 [https://www.illumina.com.cn/techniques/sequencing/rna-sequencing.html]
• 93. 生命学院杨雪瑞课题组开发单细胞剪接谱分析工具 [https://www.tsinghua.edu.cn/info/1175/110241.htm]
• 94. salmon快速基因和转录本定量| isoform [https://www.cnblogs.com/leezx/p/17673153.html]
• 95. 【亲测免费】 RSEM：RNA-Seq数据分析的强大工具原创 [https://blog.csdn.net/gitblog_01029/article/details/141237560]
• 96. 关于发布赋能药物创新的RNA基础研究重大研究计划2025 ... [https://research.bjmu.edu.cn/kyxm/gzrjjxm/kyxm_gzrjjxm_tzgg/e57168cd68b34de3ab6eb3e12f77074b.htm]
• 97. 赵方庆团队撰写人工智能驱动的环形RNA发现特邀综述 [http://www.ioz.ac.cn/gb2018/xwdt/kyjz/202504/t20250418_7603679.html]
• 98. bulk-RNA seq测序数据分析流程_bulk转录组分析 [https://blog.csdn.net/zengwanqin/article/details/135402273]
• 99. Human splicing diversity and the extent of unannotated ... [https://pmc.ncbi.nlm.nih.gov/articles/PMC5203714/]
• 100. Robust and annotation-free analysis of alternative splicing ... [https://elifesciences.org/articles/73520]
• 101. Mining alternative splicing patterns in scRNA-seq data using ... [https://genomebiology.biomedcentral.com/articles/10.1186/s13059-024-03323-6]
• 102. The utility of ultra-deep RNA sequencing in Mendelian ... [https://www.sciencedirect.com/science/article/pii/S0002929725003696]
• 103. Top Ten Ways to Improve Your RNA Isolation [https://www.thermofisher.com/us/en/home/references/ambion-tech-support/rna-isolation/general-articles/ten-ways-to-improve-your-rna-isolation.html]
• 104. Tech Tip: Tips & Best Practices for Working with RNA [https://biotium.com/tech-tips-protocols/tips-best-practices-for-working-with-rna/]
• 105. JAX_DPC: Library preparation using Watchmaker RNA ... [https://www.protocols.io/view/jax-dpc-library-preparation-using-watchmaker-rna-l-hgdnb3s5f.html]
• 106. Commentary: a review of technical considerations for planning ... [https://bmcgenomics.biomedcentral.com/articles/10.1186/s12864-025-12094-8]
• 107. Optimal selection of read depth and length for bulk RNA-Seq [https://www.revvity.com/blog/optimal-selection-read-depth-and-length-bulk-rna-seq-updated-perspective]
• 108. Detecting differential alternative events in scRNA- splicing with or... seq [https://journals.plos.org/ploscompbiol/article?id=10.1371/journal.pcbi.1007925]
• 109. How can I create a Metadata Table for my RNA-Seq ... [https://qiagen.my.salesforce-sites.com/KnowledgeBase/KnowledgeNavigatorPage?id=kA41i000000g4avCAA]
• 110. Development and evaluation of a capillary electrophoresis ... [https://pubmed.ncbi.nlm.nih.gov/40456449/]
• 111. Comparison of multiplex PCR capillary electrophoresis ... [https://www.frontiersin.org/journals/public-health/articles/10.3389/fpubh.2024.1421774/full?]
• 112. Comparisons between capillary zone electrophoresis and ... [https://pubmed.ncbi.nlm.nih.gov/16725393/]
• 113. Comparisons between capillary zone electrophoresis and ... [https://www.sciencedirect.com/science/article/abs/pii/S1570023206003515]
• 114. Evaluation of relative quantification of alternatively spliced ... [https://pmc.ncbi.nlm.nih.gov/articles/PMC5634819/]
• 115. IsoPrimer: a pipeline for designing isoform-aware primer ... [https://pmc.ncbi.nlm.nih.gov/articles/PMC12311343/]
• 116. Digital Droplet PCR for the Absolute Quantification of Exon ... [https://journals.plos.org/plosone/article?id=10.1371/journal.pone.0162467]
• 117. Identification of Splice Variants and Isoforms in ... [https://pmc.ncbi.nlm.nih.gov/articles/PMC10840079/]
• 118. Absolute Quantification of mRNA Isoforms in Adult Stem Cells ... [https://pmc.ncbi.nlm.nih.gov/articles/PMC10501916/]
• 119. High-dimensional droplet digital PCR of multiple hepatitis ... [https://pure.johnshopkins.edu/en/publications/high-dimensional-droplet-digital-pcr-of-multiple-hepatitis-b-spli/]
• 120. Droplet Digital PCR Improves Detection of BRCA1/2 Copy ... [https://pmc.ncbi.nlm.nih.gov/articles/PMC12295428/]
• 121. A systematic benchmark of Nanopore long-read RNA ... [https://www.nature.com/articles/s41592-025-02623-4]
• 122. Opportunities and challenges in long-read sequencing data ... [https://genomebiology.biomedcentral.com/articles/10.1186/s13059-020-1935-5]
• 123. Integrating short-read and long-read single-cell RNA ... [https://pubmed.ncbi.nlm.nih.gov/40050124/]
• 124. RNA Sequencing and Analysis - PMC - PubMed Central - NIH [https://pmc.ncbi.nlm.nih.gov/articles/PMC4863231/]
• 125. A comprehensive benchmarking of differential splicing ... [https://pubmed.ncbi.nlm.nih.gov/37020334/]
• 126. Alternative splicing analysis benchmark with DICAST [https://www.rna-seqblog.com/alternative-splicing-analysis-benchmark-with-dicast/]
• 127. 2025 search tool benchmark: key metrics to evaluate ... [https://www.glean.com/blog/2025-search-tool-benchmark-key-metrics-to-evaluate-accuracy-and-speed]
• 128. Benchmarking AI Agents in 2025: Top Tools, Metrics & ... [https://metadesignsolutions.com/benchmarking-ai-agents-in-2025-top-tools-metrics-performance-testing-strategies/]
• 129. Incorporating alternative splicing and mRNA editing into ... [https://pmc.ncbi.nlm.nih.gov/articles/PMC4280019/]
• 130. Integrating whole genome and transcriptome sequencing ... [https://www.nature.com/articles/s41467-025-64336-8]
• 131. Splicify – a proteogenomic pipeline to detect differentially ... [https://www.rna-seqblog.com/splicify-a-proteogenomic-pipeline-to-detect-differentially-expressed-protein-isoforms-based-on-integrating-rna-seq-and-mass-spec-data/]
• 132. Bulk RNAseq Methodology and Analysis Guide [https://www.cores.emory.edu/eicc/_includes/documents/sections/resources/RNAseq_Methodology.html]
• 133. Introduction to Bulk RNAseq data analysis [https://bioinformatics-core-shared-training.github.io/cruk-summer-school-2021/RNAseq/Markdowns/07_Data_Exploration.html]