Overview

SpliceImpactR links differential splicing to protein-level outcomes, including sequence changes, domain gain/loss, and potential PPI rewiring.

SpliceImpactR is a self-contained and efficient pipeline designed to identify how alternative RNA processing events affect resulting proteins. After initial alternative RNA processing quantification, the method runs entirely within R and does not require additional external software.

The pipeline begins by importing alternative RNA processing events that differ across conditions. It is built to support highly flexible input, requiring only inclusion and exclusion coordinates for each event isoform. This makes it compatible with both user-defined events and outputs from existing tools.

Custom input functions detect multiple classes of alternative RNA processing, including alternative first exons (AFE), hybrid first exons (HFE), retained introns (RI), skipped exons (SE), alternative 5′ and 3′ splice sites (A5SS and A3SS), mutually exclusive exons (MXE), hybrid last exons (HLE), and alternative last exons (ALE).

SpliceImpactR then performs differential inclusion analysis to identify condition-specific transcript pairs that differ by at least one alternative RNA processing event. These events are mapped to transcript and protein annotations to evaluate their effects on coding potential, sequence similarity, and frameshifts.

The pipeline next uses curated and user-provided protein feature sets to annotate domain-level changes, assess shifts in domain usage across conditions, and integrate domain-domain interaction data, domain-motif interaction data, and protein-protein interaction information to predict downstream effects on protein interaction networks.

Finally, SpliceImpactR performs a global analysis to detect co-regulated splicing events and quantify the relative contribution of each alternative RNA processing type. The pipeline is also available through an R Shiny interface to improve accessibility and ease of use.

SpliceImpactR is built to handle both S4 SpliceImpactResult object with custom accession functions and standard data.table input / workflow. SpliceImpactResult objects contain all relevant data along with detailed sequence information for each event. This object allows for maximal interoperability with any level of data from this pipeline.

This vignette uses bundled test data so every core step is runnable.

Installation

Install with one of the following methods.

Option 1: BiocManager

if (!requireNamespace("BiocManager", quietly = TRUE)) {
  install.packages("BiocManager")
}
BiocManager::install("SpliceImpactR")

Option 2: devtools (GitHub)

if (!requireNamespace("devtools", quietly = TRUE)) {
  install.packages("devtools")
}
devtools::install_github("fiszbein-lab/SpliceImpactR")

External tools and acronyms

rMATS (replicate Multivariate Analysis of Transcript Splicing) provides event-level splicing tables.
HITindex quantifies first/last exon usage and labels hybrid exons.
PPIDM Protein-Protein Interaction Domain Miner for domain-domain Interaction derived from PPI and 3did’s domain-domain interactions
ELM Eukaryotic Linear Motif Database for accessing Short Linear Motif occurences and domain-motif interactions
BiomaRt accesses data from InterPro, PFAM, SignalP, TMHMM, CDD, Mobidb-lite

Workflow map

Load reference resources (annotations + protein features).
Build sample manifest (sample_frame).
(Optional) run the quick-start wrapper for end-to-end execution.
Read splicing events for stepwise analysis.
Run QC summaries and plots.
Run differential inclusion.
Select case/control transcript pairs for DI events (get_ranked_pairs()).
Quantify sequence and domain impacts.
Infer PPI rewiring.
Build enrichment foreground gene sets from table or S4 inputs.
Visualize transcript and protein consequences for one event pair.
Generate integrated multi-layer summaries.
Use S4 container workflows and accessors.
Run optional custom-input entry points.

Load Reference Resources

get_protein_features() supports these feature sources: - interpro: integrated domain/family/superfamily signatures. - pfam: protein domain HMM families. - cdd: NCBI Conserved Domain Database annotations. - gene3d: CATH/Gene3D structural domain annotations. - signalp: signal peptide calls. - tmhmm: transmembrane helix calls. - ncoils: coiled-coil region predictions. - seg: low-complexity segments. - mobidblite: intrinsically disordered region calls. - elm: short linear motifs (SLiMs; retrieved from ELM and mapped to transcripts through BioMart).

For BioMart-backed sources, each feature must provide three attributes: {feature}, {feature}_start, and {feature}_end (for example pfam, pfam_start, pfam_end). You can add non-BioMart custom features with get_manual_features() (explained further down in the Custom Input Entry Points section) and merge everything through get_comprehensive_annotations().

Full runs query Ensembl BioMart for release (default 111, the Ensembl release of the default GENCODE v45 annotation); pass the Ensembl release that matches your get_annotation() release. Ensembl 116 (June 2026) is the last release with BioMart, so each release is queried at its own address (for example https://e111.ensembl.org, which Ensembl forwards to the 111 archive), and ensembl_host selects another BioMart host.

We load test annotations and three representative feature sources (interpro, signalp, elm) and combine them. Loading multiple databases in one call is possible, but not recommended due to timeouts. These are automatically cached and loaded from cache if the correct signature (organism, version, etc) is found. Ensembl BioMart is sometimes temporarily unavailable (HTTP 500/503 or time-outs), so this procedure will sometimes take a few tries. If this is unsuccessful, try again later, choose another release, restart R and/or try suggested methods related to this issue (eg: httr::reset_config())

library(SpliceImpactR)
annotation_df <- get_annotation(load = "test")
#> [INFO] Loading bundled test annotation data

interpro_features <- get_protein_features(
  biomaRt_databases = c("interpro"),
  gtf_df = annotation_df$annotations,
  test = TRUE
)

signalp_features <- get_protein_features(
  biomaRt_databases = c("signalp"),
  gtf_df = annotation_df$annotations,
  test = TRUE
)

elm_features <- get_protein_features(
  biomaRt_databases = c("elm"),
  gtf_df = annotation_df$annotations,
  test = TRUE
)

protein_feature_total <- get_comprehensive_annotations(
  list(interpro_features, signalp_features, elm_features)
)

exon_features <- get_exon_features(
  annotation_df$annotations,
  protein_feature_total
)

For a full reference, get_annotation(load = "link") and load = "path" default to TSL 1, 2, and 3. Set filter_tsl = NULL to retain transcripts regardless of their support level, including missing or unknown TSL values. Use filter_tsl = NULL again with load = "cached" to retrieve that reference; it is stored separately from TSL-filtered references. Gene-biotype and incomplete-CDS tag filters still apply. The bundled test reference is loaded as stored.

table(protein_feature_total$database)
#> 
#>      elm interpro  signalp 
#>       18     1449       13

Build Sample Manifest

Each sample directory contains event files exported by rMATS and/or HITindex. The manifest below points to the bundled test directories. Using this function, SpliceImpactR will search each directory for the outputs provided by rMATS and HIT Index output files ending in .AFEPSI, .ALEPSI, and .exon. Label one condition case and one condition control.

Required sample manifest columns: - path: per-sample directory containing splicing output files. - sample_name: unique sample identifier. - condition: case/control label used in DI and downstream pairing.

Manifest expectations: - One row per sample. - path values must point to readable sample directories. - Use exactly two conditions for this standard workflow (case, control). - Include replicate samples per condition for robust differential modeling.

sample_frame <- data.frame(
  path = c(
    file.path(system.file("extdata", package = "SpliceImpactR"),
              "rawData/control_S5/"),
    file.path(system.file("extdata", package = "SpliceImpactR"),
              "rawData/control_S6/"),
    file.path(system.file("extdata", package = "SpliceImpactR"),
              "rawData/control_S7/"),
    file.path(system.file("extdata", package = "SpliceImpactR"),
              "rawData/control_S8/"),
    file.path(system.file("extdata", package = "SpliceImpactR"),
              "rawData/case_S1/"),
    file.path(system.file("extdata", package = "SpliceImpactR"),
              "rawData/case_S2/"),
    file.path(system.file("extdata", package = "SpliceImpactR"),
              "rawData/case_S3/"),
    file.path(system.file("extdata", package = "SpliceImpactR"),
              "rawData/case_S4/")
  ),
  sample_name = c("S5", "S6", "S7", "S8", "S1", "S2", "S3", "S4"),
  condition = c(
    "control", "control", "control", "control",
    "case", "case", "case", "case"
  ),
  stringsAsFactors = FALSE
)

Quick start wrapper (`get_splicing_impact`)

If your inputs already follow the standard layout, run the full workflow with the top-level wrapper and choose table or S4 return mode.

# data.table-first return
out_dt <- get_splicing_impact(
  sample_frame = sample_frame,
  source_data = "both",
  annotation_df = annotation_df,
  protein_feature_total = protein_feature_total,
  exon_features = exon_features,
  return_class = "data.table"
)

# S4-first return (single object carries all slots)
out_s4 <- get_splicing_impact(
  sample_frame = sample_frame,
  source_data = "both",
  annotation_df = annotation_df,
  protein_feature_total = protein_feature_total,
  exon_features = exon_features,
  return_class = "S4"
)

The wrapper selects transcript pairs with the ORF-aware matcher by default (matching = "orf", described under transcript pairs below); matching = "legacy" runs the previous matcher. With the ORF matcher, out_dt$matching holds the candidate, ranking and unresolved-comparison diagnostics.

Read Splicing Events

For the stepwise workflow in this vignette, load event data explicitly:

data <- get_rmats_hit(
  sample_frame,
  event_types = c("ALE", "AFE", "MXE", "SE", "A3SS", "A5SS", "RI")
)
#> [INFO] Filtered 0 unannotated ALE/AFE events (465 -> 465)

DT <- data.table::as.data.table(data)
DT[, .(
  n_rows = .N,
  n_events = data.table::uniqueN(event_id),
  n_genes = data.table::uniqueN(gene_id)
)]
#>    n_rows n_events n_genes
#>     <int>    <int>   <int>
#> 1:   5905      614      41

If rMATS and HITindex outputs are stored in separate per-sample directories, load them independently and combine on shared columns. In this test-data example, both manifests point to the same bundled paths; in real analyses, sample_frame_rmats$path and sample_frame_hit$path can be different.

sample_frame_rmats <- sample_frame
sample_frame_hit <- sample_frame

rmats_only <- get_rmats(
  load_rmats(
    sample_frame_rmats,
    use = "JCEC",
    event_types = c("MXE", "SE", "A3SS", "A5SS", "RI")
  )
)

hit_only <- get_hitindex(
  sample_frame_hit,
  keep_annotated_first_last = TRUE
)
#> [INFO] Filtered 0 unannotated ALE/AFE events (465 -> 465)

shared_cols <- intersect(names(rmats_only), names(hit_only))
data_separate <- data.table::rbindlist(
  list(
    rmats_only[, ..shared_cols],
    hit_only[, ..shared_cols]
  ),
  use.names = TRUE,
  fill = TRUE
)

data_separate[, .(
  n_rows = .N,
  n_events = data.table::uniqueN(event_id),
  n_genes = data.table::uniqueN(gene_id)
)]
#>    n_rows n_events n_genes
#>     <int>    <int>   <int>
#> 1:   5905      614      41

Sample-Level QC and Exploration

Compare HITindex Between Conditions

This summary compares event-level mean HIT values across conditions.

hit_compare <- compare_hit_index(
  sample_frame,
  condition_map = c(control = "control", test = "case")
)

hit_compare$plot

HITindex condition comparison.

head(hit_compare$results[order(fdr)], 6)
#>                               event_key    control        test   diff_HIT
#>                                  <char>      <num>       <num>      <num>
#> 1: ENSG00000158286|chr1:6210370-6211118 -0.1890000  0.12766667  0.3166667
#> 2: ENSG00000158286|chr1:6211867-6212053 -0.1905000 -0.23400000 -0.0435000
#> 3: ENSG00000158286|chr1:6212218-6212416 -0.1783333 -0.14533333  0.0330000
#> 4: ENSG00000158286|chr1:6212682-6213183  0.4350000  0.04666667 -0.3883333
#> 5: ENSG00000158286|chr1:6218289-6218369 -0.2993333 -0.05766667  0.2416667
#> 6: ENSG00000142599|chr1:8355419-8355599  0.1010000  0.10800000  0.0070000
#>    delta_HIT n_control n_test   p_value   fdr
#>        <num>     <int>  <int>     <num> <num>
#> 1: 0.3166667         2      3 1.0000000     1
#> 2: 0.0435000         2      3 1.0000000     1
#> 3: 0.0330000         3      3 0.8337307     1
#> 4: 0.3883333         2      3 1.0000000     1
#> 5: 0.2416667         3      3 0.2567179     1
#> 6: 0.0070000         4      4 0.9338296     1

PSI Distribution Overview

This panel summarizes event count depth-normalized metrics and PSI eCDF.

overview_plot <- overview_spicing_comparison(
  events = data,
  sample_df = sample_frame,
  depth_norm = "exon_files",
  event_type = "AFE"
)

overview_plot

Depth-normalized event counts and PSI overview.

Differential Inclusion (DI)

We run quasi-binomial GLMs and filter significant events. Verbose filtering messages are line-wrapped so progress logs remain readable without horizontal scrolling. For efficient processing, set parallel_glm to TRUE and provide a BiocParallel object to BPPARAM.

Default thresholds in keep_sig_pairs is FDR < 0.05 and |DELTA PSI| > 0.1.

res <- get_differential_inclusion(
  data,
  min_total_reads = 10,
  parallel_glm = TRUE,
  BPPARAM = BiocParallel::SerialParam()
)
#> [INFO] Input contains 41 genes, 614 events, and 2952 event instances.
#> [PROCESSING/INFO] Filtering out low-coverage rows removed 914 event instances; 4991
#> remain.
#> [PROCESSING/INFO] Completing AFE/ALE with zeros per sample (total strategy: event_max)
#> filled 407 AFE rows and 160 ALE rows (totals: AFE=608, ALE=424).
#> [PROCESSING/INFO] Filtering genes with no nonzero values per sample/event_type removed
#> 37 genes from specific sample groups; 40 genes remain overall.
#> [PROCESSING/INFO] Filtering by minimum condition presence removed 213 events; 354
#> events remain.
#> [PROCESSING/INFO] Removed 63 non-changing events; 291 events remain.
#> [PROCESSING/INFO] Filtering events not represented in both conditions removed 32
#> events; 259 events remain.
#> [PROCESSING] Fitting quasi-binomial GLMs per site...
#>   |                                                                              |                                                                      |   0%  |                                                                              |======================================================================| 100%
#> 
#> [DONE] Fitted 532 sites in 1 chunks.
#> [INFO] Done.

res_di <- keep_sig_pairs(res)

res[, .(
  n_rows = .N,
  n_sig = sum(padj <= 0.05 & abs(delta_psi) >= 0.1, na.rm = TRUE)
)]
#>    n_rows n_sig
#>     <int> <int>
#> 1:    533    72

plot_di_volcano_dt(res)

Differential inclusion volcano plot.

Match Events to Sequences and Build Transcript Pairs

For each significant comparison, get_ranked_pairs() selects one case and one control transcript that represent the two forms and differ as little as possible outside the event (the default matcher in get_splicing_impact()):

Validate the input: missing metadata, strands other than +/-, unknown forms and malformed intervals stop the run, naming the affected events.
Pair forms: positive delta_psi rows are the case and negative rows the control, for each event_id.
Keep, for each form, same-gene transcripts that use the event’s defining splice sites exactly and have no exon in the excluded region. Outer flank ends and TSS/PAS ends may differ (labelled relaxed). A form with more than max_candidates (100) candidates keeps the best supported, and the number dropped is reported.
Rank candidate pairs by coding status (an annotated CDS, so NMD-annotated transcripts count), transcript support tier (TSL 1-3, 4-5, unknown), how much coding-region structure the two transcripts share outside the event (a frameshifted isoform is not penalised for its early stop), and then exact TSL. A transcript is never paired with itself.
When no structural pair exists, choose a pair with the legacy overlap rules and label it approximate (fallback = FALSE disables this). Comparisons still without a pair are listed in unmatched with a reason.

For an event with several sites, such as an AFE with three first exons, each rising x falling pair of sites is its own comparison, so a site can appear in several rows. keep_sig_pairs() keeps every site of an event in which any site passes; site_significant_case/_control show which sites passed themselves, and n_event_comparisons gives the number of comparisons the event defines. Pair counts, including domain enrichment, count comparisons rather than events.

ranked <- get_ranked_pairs(
  res_di,
  annotation_df$annotations,
  annotation_df$sequences
)
#> [RANK] Form pair 1/52
#> [RANK] Form pair 25/52
#> [RANK] Form pair 50/52
pairs <- ranked$pairs
hits_sequences <- ranked$matched

pairs[, .(
  n_pairs = .N,
  n_events = data.table::uniqueN(event_id)
)]
#>    n_pairs n_events
#>      <int>    <int>
#> 1:      33       26
table(pairs$matching_tier)
#> 
#> approximate  structural 
#>          17          16

With the small test annotation, many comparisons lack a transcript with the exact event structure, so about half the pairs are approximate; with a full GENCODE annotation most pairs are structural.

ranked also holds every candidate with its eligibility, every scored pair, the unresolved comparisons and the settings. The transcript_pair_ranking vignette describes the protocol step by step; scores and ties are descriptive diagnostics, not expression probabilities.

Legacy matcher

The previous matcher (matching = "legacy" in get_splicing_impact()) picks each form’s transcript independently. It is faster, but it does not check the event’s exact structure and can pair a transcript with itself. It works through a strict hierarchy:

Prefilter by chr, strand, and gene_id to keep only compatible annotation intervals.
Match inclusion (inc) coordinates to exons and require all inclusion parts to be covered for a transcript candidate.
Remove candidates that overlap exclusion (exc) coordinates.
Prioritize transcript/exon choices by event-type-consistent exon class (first, internal, last), then reciprocal overlap and intersection width; protein-coding status and then TSL only break remaining ties.
Build case/control pairs in get_pairs(source = "multi") by joining all positive delta_psi rows (case) with all negative rows (control) for each event_id, then ordering by strongest |delta_psi|. Forms that match no transcript are left out.

matched <- get_matched_events_chunked(
  res_di,
  annotation_df$annotations,
  chunk_size = 2000
)
hits_sequences <- attach_sequences(matched, annotation_df$sequences)
pairs <- get_pairs(hits_sequences, source = "multi")

Once pairing is done, we can probe whether there is a proximal/distal shift in AFE/ALE usage.

proximal_output <- get_proximal_shift_from_hits(pairs)
proximal_output$plot

Proximal versus distal terminal exon summary.

Sequence-Level Consequences

compare_sequence_frame() computes protein-coding class, alignment identity, frameshift/rescue status, and length deltas.

Key labels used downstream: - protein_coding: both isoforms have protein IDs (both sides coding). - onePC: only one side has a protein ID. - noPC: neither side has a protein ID. - Match: protein sequences are identical. - FrameShift: frame is disrupted between paired isoforms. - NMD: an isoform is annotated as nonsense-mediated decay, so it is taken to make no protein (frame_call still shows the frameshift that usually causes it).

seq_compare <- compare_sequence_frame(pairs, annotation_df$annotations)
#> [1] "[INFO] Processing 33 transcript and protein sequence alignments, this may take a little..."
#> [1] "[Processing] Identifying frame shifts and rescues"
#> [1] "[INFO] 1 frameshifts (0 rescues) and 16 non-frameshifts were identified, "

alignment_plot <- plot_alignment_summary(seq_compare)
alignment_plot

Protein alignment summary by event type.

length_plot <- plot_length_comparison(seq_compare)
length_plot

Case versus control isoform length comparison.

Domain-Level Changes and Domain Enrichment

Background transcript pairs are derived from annotations and compared against observed event-linked domain changes. source = "annotated" is the default. A HIT-index background (source = "hit_index") keeps one transcript per observed exon, so it usually lacks the matched pairs and domain enrichment rejects it; its gene universes still work for gene-level enrichment.

bg <- get_background(
  source = "annotated",
  annotations = annotation_df$annotations,
  protein_features = protein_feature_total,
  BPPARAM = BiocParallel::SerialParam()
)

hits_domain <- get_domains(seq_compare, exon_features)

hits_domain[, .(
  n_hits = .N,
  n_domain_changed = sum(diff_n > 0, na.rm = TRUE)
)]
#>    n_hits n_domain_changed
#>     <int>            <int>
#> 1:     33               10

enriched_domains <- enrich_domains_hypergeo(
  hits = hits_domain,
  background = bg,
  db_filter = "interpro"
)

domain_plot <- plot_enriched_domains_counts(enriched_domains, top_n = 20)
domain_plot

Top enriched InterPro domain differences.

You can also stratify enrichment by event type or by feature database.

enriched_afe_interpro <- enrich_by_event(
  hits = hits_domain,
  background = bg,
  events = "AFE",
  db_filter = "interpro"
)

enriched_interpro_only <- enrich_by_db(
  hits = hits_domain,
  background = bg,
  dbs = "interpro"
)

list(
  by_event_rows = nrow(enriched_afe_interpro),
  by_db_rows = nrow(enriched_interpro_only)
)
#> $by_event_rows
#> [1] 1
#> 
#> $by_db_rows
#> [1] 1

PPI Rewiring

PPI effects are inferred from domain and motif changes. Domain-domain and domain-motif interactions are sourced from 3did/PPDIM and ELM. PPI rewiring currently only works for human inputs.

ppi <- get_ppi_interactions()

hits_final <- get_ppi_switches(
  hits_domain,
  ppi,
  protein_feature_total
)
#> 

ppi_plot <- plot_ppi_summary(hits_final)
ppi_plot

Summary of predicted case/control PPI changes.

Gene Enrichment Foreground Probes

get_gene_enrichment() extracts foreground genes for enrichment.

mode = "di" expects a DI-like table (res).
mode = "ppi"/"domain" expects a hits-like table (hits_final or hits_domain).
S4 examples are shown in Section 13 so object setup stays in one place.

fg_di <- get_gene_enrichment(
  mode = "di",
  res = res,
  padj_threshold = 0.1,
  delta_psi_threshold = 0.1
)

fg_domain <- get_gene_enrichment(mode = "domain", hits = hits_domain)
fg_ppi <- get_gene_enrichment(mode = "ppi", hits = hits_final)

lengths(list(di = fg_di, domain = fg_domain, ppi = fg_ppi))
#>     di domain    ppi 
#>     39      9      5
intersect(fg_di, fg_ppi)
#> [1] "ENSG00000142599" "ENSG00000140983" "ENSG00000101391" "ENSG00000116406"

The optional get_enrichment() step depends on additional packages and may return no significant terms on small toy datasets. If you see all NA statistics, foreground genes often map too sparsely into the selected background and ontology. p_adjust_cutoff and top_n_plot are set to allow for example visualization, but should be adjusted to standard (eg: 0.05, 10) for non-toy data uses.

needed <- c("clusterProfiler", "msigdbr", "org.Hs.eg.db")
has_needed <- all(vapply(
  needed,
  requireNamespace,
  FUN.VALUE = logical(1),
  quietly = TRUE
))
#> 
#> 

if (has_needed) {
  enrichment_di <- get_enrichment(
    foreground = fg_domain,
    background = intersect(attr(bg, "feature_gene_universe"),
                           res$gene_id[is.finite(res$padj)]),
    species = "human",
    gene_id_type = "ensembl",
    sources = "GO:BP",
    p_adjust_cutoff = 1, ## CHANGE FOR NON-TOY DATA USAGE, eg: 0.05
    top_n_plot = 1 ## CHANGE FOR NON-TOY DATA USAGE, eg: NULL or 10
  )

  enrichment_di$plot
} else {
  message(
    "Skipping get_enrichment(): install ",
    paste(needed, collapse = ", ")
  )
}
#> 'select()' returned 1:1 mapping between keys and columns
#> 'select()' returned 1:1 mapping between keys and columns
#> Warning: semData is not provided. It will be calculated automatically.

Visualize Transcript and Protein Consequences

plot_two_transcripts_with_domains_unified() shows matched isoforms in either transcript (genome-coordinate) or protein (amino-acid) view.

viz_pair <- hits_final[
  !is.na(transcript_id_case) &
    !is.na(transcript_id_control) &
    transcript_id_case != "" &
    transcript_id_control != ""
][1]

if (nrow(viz_pair) == 0) {
  stop("No transcript pairs available for visualization.")
}

tx_pair <- c(viz_pair$transcript_id_case, viz_pair$transcript_id_control)

transcript_centric <- plot_two_transcripts_with_domains_unified(
  transcripts = tx_pair,
  gtf_df = annotation_df$annotations,
  protein_features = protein_feature_total,
  feature_db = c("interpro"),
  combine_domains = TRUE,
  view = "transcript"
)

protein_centric <- plot_two_transcripts_with_domains_unified(
  transcripts = tx_pair,
  gtf_df = annotation_df$annotations,
  protein_features = protein_feature_total,
  feature_db = c("interpro"),
  combine_domains = TRUE,
  view = "protein"
)

transcript_centric

Transcript and protein views for one paired event.

protein_centric

Transcript and protein views for one paired event.

Inspect a Single Event

probe_individual_event() plots PSI for one event across samples.

event_to_probe <- res$event_id[1]
probe <- probe_individual_event(data, event = event_to_probe)
probe$plot

PSI distribution for one selected event.

head(probe$data)
#>    event_id sample condition   psi i.condition   form
#>      <char> <char>    <char> <num>      <char> <char>
#> 1:  A3SS:28     S1      case 0.483        case    INC
#> 2:  A3SS:28     S2      case 0.483        case    INC
#> 3:  A3SS:28     S3      case 0.525        case    INC
#> 4:  A3SS:28     S4      case 0.345        case    INC
#> 5:  A3SS:28     S5   control 0.208     control    INC
#> 6:  A3SS:28     S6   control 0.085     control    INC

Integrative Visualization

integrated_event_summary() combines sequence, domain, and PPI summaries.

Panel guide: - Top-left: event classification composition (Match, FrameShift, etc.) by event type. - Top-right: alignment score distributions by event type. - Middle-left: domain-change prevalence (Any, case-only, control-only, both). - Middle-center/right: fraction and counts of predicted PPI rewiring. - Bottom-left: relative event retention from pre-filter DI to final hits. - Bottom-right: gene-level event-type coordination (Jaccard heatmap).

int_summary <- integrated_event_summary(hits_final, res)
int_summary$plot

Integrated multi-layer summary.

int_summary$summaries$relative_use
#>    event_type n_pre n_post relative_use
#>        <fctr> <int>  <int>        <num>
#> 1:        AFE    22     15   0.68181818
#> 2:       A3SS    34      5   0.14705882
#> 3:         RI    35      2   0.05714286
#> 4:         SE   135      4   0.02962963
#> 5:       A5SS    12      0   0.00000000
#> 6:        ALE    10      0   0.00000000
#> 7:        MXE    11      0   0.00000000

Output Description (`hits_final`, `data`, `res`)

The pipeline returns three core tables in data.table mode: - data: sample-level event measurements before differential modeling. - res: differential inclusion results per tested event/site. - hits_final: paired case/control isoform effects with sequence, domain, and PPI annotations.

Suffix convention used throughout: - _case = case-preferred isoform values. - _control = control-preferred isoform values.

`hits_final` (integrated event-level output)

Use this table for biological interpretation and downstream plotting.

Event and isoform identifiers: event_id, event_type, gene_id, chr, strand, transcript_id_case, transcript_id_control, protein_id_case, protein_id_control, form_case, form_control, exons_case, exons_control.
Event coordinates and DI statistics: inc_case, inc_control, exc_case, exc_control, delta_psi_case, delta_psi_control, p.value_case, p.value_control, padj_case, padj_control, n_samples_case, n_samples_control, n_case_case, n_case_control, n_control_case, n_control_control.
Sequence and coding context: transcript_seq_case, transcript_seq_control, protein_seq_case, protein_seq_control, pc_class, length columns (prot_len_*, tx_len_*, exon_cds_len_*, exon_len_*) and their *_diff / *_diff_abs derivatives.
Alignment/frame classification: DNA alignment (dna_pid, dna_score, dna_width), protein alignment (prot_pid, prot_score, prot_width), frame diagnostics (frame_call, rescue, frame_check_exon1, frame_check_exon2), and summary_classification.
Domain-level change annotations: domains_exons_case, domains_exons_control, case_only_domains, control_only_domains, case_only_domains_list, control_only_domains_list, either_domains_list, case_only_n, control_only_n, diff_n.
Predicted interaction rewiring: case_ppi, control_ppi, n_case_ppi, n_control_ppi, n_ppi, case_ppi_drivers, control_ppi_drivers.

`data` (raw sample-level input table)

Use data to inspect per-sample evidence feeding DI.

Core columns: event_id, event_type, form, gene_id, chr, strand, inc, exc, inclusion_reads, exclusion_reads, psi, sample, condition, source_file.

Often present depending on import path: HITindex metadata such as HITindex, class, nFE, nLE, nUP, nDOWN, nTXPT, psi_original, total_reads, source.

`res` (differential inclusion output)

Use res to rank significant events before pairing/domain/PPI steps.

Core columns: site_id, event_id, event_type, gene_id, chr, strand, inc, exc, form, n_samples, n_control, n_case, mean_psi_ctrl, mean_psi_case, delta_psi, p.value, padj, cooks_max, n, n_used.

S4 Applications and Accessors

Use SpliceImpactResult when you want all pipeline pieces in one object, with consistent filtering across event-linked slots.

SpliceImpactResult is a custom S4 container that keeps all major pipeline parts synchronized: - raw_events (SummarizedExperiment): sample-level table + ranges/assays. - di_events / res_di (GRanges): differential inclusion rows. - matched (DFrame): annotation-matched rows (and sequence-attached rows). - paired_hits (GRanges) + segments (GRangesList): final case/control event impacts. - sample_frame (DFrame): sample manifest metadata.

obj <- as_splice_impact_result(
  data = data,
  res = res,
  res_di = res_di,
  matched = hits_sequences,
  sample_frame = sample_frame,
  hits_final = hits_final
)

get_hits_*() wrappers provide compact views for common downstream tasks: - get_hits_core(): event IDs, transcript/protein IDs, and key DI metadata. - get_hits_domain(): domain gain/loss columns (case_only_*, control_only_*, diff_n, mapped domain lists). - get_hits_ppi(): PPI partner changes (n_ppi, case/control partner lists, PPI drivers). - get_hits_sequence(): sequence/alignment/frame fields (pc_class, frame_call, protein/transcript length deltas, identity metrics).

raw_dt <- as_dt_from_s4(obj, "raw_events")
di_dt <- as_dt_from_s4(obj, "di_events")
hits_dt <- as_dt_from_s4(obj, "paired_hits")

hits_core <- get_hits_core(obj)
hits_domain_view <- get_hits_domain(obj)
hits_ppi_view <- get_hits_ppi(obj)
hits_sequence_view <- get_hits_sequence(obj)

c(
  raw_rows = nrow(raw_dt),
  di_rows = nrow(di_dt),
  hits_rows = nrow(hits_dt),
  hits_core_rows = nrow(hits_core),
  hits_domain_rows = nrow(hits_domain_view),
  hits_ppi_rows = nrow(hits_ppi_view),
  hits_sequence_rows = nrow(hits_sequence_view)
)
#>           raw_rows            di_rows          hits_rows     hits_core_rows 
#>               5905                533                 33                 33 
#>   hits_domain_rows      hits_ppi_rows hits_sequence_rows 
#>                 33                 33                 33

obj_focus <- filter_spliceimpact_hits(
  obj,
  diff_n > 0,
  n_ppi > 0
)

focus_hits <- as_dt_from_s4(obj_focus, "paired_hits")
focus_hits[, .N, by = event_type][order(-N)]
#>    event_type     N
#>        <char> <int>
#> 1:         RI     2
#> 2:         SE     2
#> 3:        AFE     1

fg_di_s4 <- get_gene_enrichment(
  mode = "di",
  x = obj,
  padj_threshold = 0.05,
  delta_psi_threshold = 0.1
)

fg_ppi_s4 <- get_gene_enrichment(mode = "ppi", x = obj)

fg_ppi_focus <- get_gene_enrichment(mode = "ppi", x = obj_focus)
list(
  di_equal_table = identical(sort(fg_di), sort(fg_di_s4)),
  ppi_equal_table = identical(sort(fg_ppi), sort(fg_ppi_s4)),
  focused_ppi_genes = head(fg_ppi_focus)
)
#> $di_equal_table
#> [1] TRUE
#> 
#> $ppi_equal_table
#> [1] TRUE
#> 
#> $focused_ppi_genes
#> [1] "ENSG00000142599" "ENSG00000140983" "ENSG00000164692" "ENSG00000116406"
#> [5] "ENSG00000101391"

spliceimpact_s4_schema()
#> $slots
#> $slots$raw_events
#> [1] "SummarizedExperiment: sample/form rows with rowRanges + rowData"
#> 
#> $slots$di_events
#> [1] "GRanges: differential inclusion rows with mcols"
#> 
#> $slots$res_di
#> [1] "GRanges: threshold-filtered differential rows with mcols"
#> 
#> $slots$matched
#> [1] "S4Vectors::DataFrame: annotation-matched DI rows (and sequence-attached rows when present)"
#> 
#> $slots$sample_frame
#> [1] "S4Vectors::DataFrame: sample manifest (`path`, `sample_name`, `condition`)"
#> 
#> $slots$paired_hits
#> [1] "GRanges: paired case/control rows with mcols"
#> 
#> $slots$segments
#> [1] "GRangesList: per-pair genomic segments (inc_case/inc_control/exc_case/exc_control)"
#> 
#> $slots$metadata
#> [1] "list: provenance and pipeline metadata"
#> 
#> 
#> $keys
#> $keys$raw_events
#> [1] "raw_key (rowData)"
#> 
#> $keys$di_events
#> [1] "di_key (mcols)"
#> 
#> $keys$res_di
#> [1] "di_key (mcols)"
#> 
#> $keys$paired_hits
#> [1] "pair_key (mcols)"
#> 
#> $keys$segments
#> [1] "names(segments) == pair_key"
#> 
#> 
#> $assays
#> [1] "psi"             "inclusion_reads" "exclusion_reads"

S4-first Main Workflow (Code Only)

This mirrors the table workflow, but keeps everything in a single SpliceImpactResult object. get_ranked_pairs() stores the selected pairs and form rows in the object and its diagnostics in metadata$matching_diagnostics. For the legacy matcher, replace that step with get_matched_events_chunked(), attach_sequences() and get_pairs(source = "multi"), each with return_class = "S4".

obj_flow <- as_splice_impact_result(
  data = data,
  sample_frame = sample_frame
)

obj_flow <- get_differential_inclusion(
  obj_flow,
  min_total_reads = 10,
  parallel_glm = TRUE,
  BPPARAM = BiocParallel::SerialParam(),
  return_class = "S4"
)
obj_flow <- keep_sig_pairs(obj_flow, return_class = "S4")

obj_flow <- get_ranked_pairs(
  obj_flow,
  annotation_df$annotations,
  annotation_df$sequences,
  return_class = "S4"
)
obj_flow <- compare_sequence_frame(
  obj_flow,
  annotation_df$annotations,
  return_class = "S4"
)
obj_flow <- get_domains(obj_flow, exon_features, return_class = "S4")
obj_flow <- get_ppi_switches(
  obj_flow,
  ppi,
  protein_feature_total,
  return_class = "S4"
)

hits_core_flow <- get_hits_core(obj_flow)
hits_domain_flow <- get_hits_domain(obj_flow)
hits_ppi_flow <- get_hits_ppi(obj_flow)
hits_sequence_flow <- get_hits_sequence(obj_flow)

Custom Input Entry Points

These helpers support non-standard inputs at multiple stages of the workflow.

Add Manual Protein Features

get_manual_features() lets you add user-defined peptide-level features.

ann_dt <- data.table::as.data.table(annotation_df$annotations)
coding_tx <- unique(ann_dt[type == "exon" & cds_has == TRUE, transcript_id])
n_manual <- min(3L, length(coding_tx))

if (n_manual < 1L) {
  stop("No coding transcripts found for manual feature example.")
}

manual_df <- data.frame(
  ensembl_transcript_id = coding_tx[seq_len(n_manual)],
  ensembl_peptide_id = rep(NA_character_, n_manual),
  name = paste0("demo_feature_", seq_len(n_manual)),
  start = c(20L, 45L, 80L)[seq_len(n_manual)],
  stop = c(35L, 58L, 92L)[seq_len(n_manual)],
  database = rep("manual", n_manual),
  alt_name = rep(NA_character_, n_manual),
  feature_id = rep(NA_character_, n_manual),
  stringsAsFactors = FALSE
)

manual_features <- get_manual_features(
  manual_features = manual_df,
  gtf_df = annotation_df$annotations
)
#> Placed 3 of 3 manual feature row(s).

head(manual_features)
#>    ensembl_transcript_id feature_id start  stop    chr strand     clean_name
#>                   <char>     <char> <num> <num> <char> <char>         <char>
#> 1:       ENST00000337907       <NA>    20    35   chr1      - demo_feature_1
#> 2:       ENST00000400908       <NA>    45    58   chr1      - demo_feature_2
#> 3:       ENST00000476556       <NA>    80    92   chr1      - demo_feature_3
#>    alt_name database ensembl_peptide_id method
#>      <char>   <char>             <char> <char>
#> 1:     <NA>   manual               <NA> manual
#> 2:     <NA>   manual               <NA> manual
#> 3:     <NA>   manual               <NA> manual
#>                                   name
#>                                 <char>
#> 1: demo_feature_1;chr1:8656240-8656193
#> 2: demo_feature_2;chr1:8656165-8656124
#> 3: demo_feature_3;chr1:8362685-8361841

Pre-DI Custom Table

Use get_user_data() when you already have sample-level reads/PSI.

example_df <- data.frame(
  event_id = rep("A3SS:1", 8),
  event_type = "A3SS",
  form = rep(c("INC", "EXC"), each = 4),
  gene_id = "ENSG00000158286",
  chr = "chrX",
  strand = "-",
  inc = c(rep("149608626-149608834", 4), rep("149608626-149608829", 4)),
  exc = c(rep("", 4), rep("149608830-149608834", 4)),
  inclusion_reads = c(30, 32, 29, 31, 2, 3, 4, 3),
  exclusion_reads = c(1, 1, 2, 1, 28, 27, 26, 30),
  sample = c("S1", "S2", "S3", "S4", "S1", "S2", "S3", "S4"),
  condition = rep(c("case", "case", "control", "control"), 2),
  stringsAsFactors = FALSE
)
example_df$psi <- example_df$inclusion_reads /
  (example_df$inclusion_reads + example_df$exclusion_reads)

user_data <- get_user_data(example_df)
head(user_data)
#>    event_id event_type   form         gene_id    chr strand                 inc
#>      <char>     <char> <char>          <char> <char> <char>              <char>
#> 1:   A3SS:1       A3SS    INC ENSG00000158286   chrX      - 149608626-149608834
#> 2:   A3SS:1       A3SS    INC ENSG00000158286   chrX      - 149608626-149608834
#> 3:   A3SS:1       A3SS    INC ENSG00000158286   chrX      - 149608626-149608834
#> 4:   A3SS:1       A3SS    INC ENSG00000158286   chrX      - 149608626-149608834
#> 5:   A3SS:1       A3SS    EXC ENSG00000158286   chrX      - 149608626-149608829
#> 6:   A3SS:1       A3SS    EXC ENSG00000158286   chrX      - 149608626-149608829
#>                    exc inclusion_reads exclusion_reads        psi sample
#>                 <char>           <num>           <num>      <num> <char>
#> 1:                                  30               1 0.96774194     S1
#> 2:                                  32               1 0.96969697     S2
#> 3:                                  29               2 0.93548387     S3
#> 4:                                  31               1 0.96875000     S4
#> 5: 149608830-149608834               2              28 0.06666667     S1
#> 6: 149608830-149608834               3              27 0.10000000     S2
#>    condition source_file
#>       <char>      <char>
#> 1:      case            
#> 2:      case            
#> 3:   control            
#> 4:   control            
#> 5:      case            
#> 6:      case

Post-DI Custom Table

Use get_user_data_post_di() when you already have event-level DI summaries.

example_user_data <- data.frame(
  event_id = rep("A3SS:1", 8),
  event_type = "A3SS",
  gene_id = "ENSG00000158286",
  chr = "chrX",
  strand = "-",
  form = rep(c("INC", "EXC"), each = 4),
  inc = c(rep("149608626-149608834", 4), rep("149608626-149608829", 4)),
  exc = c(rep("", 4), rep("149608830-149608834", 4)),
  stringsAsFactors = FALSE
)

user_res <- get_user_data_post_di(example_user_data)
head(user_res)
#>                                                                  site_id
#>                                                                   <char>
#> 1:                    A3SS|ENSG00000158286|chrX|149608626-149608834||INC
#> 2:                    A3SS|ENSG00000158286|chrX|149608626-149608834||INC
#> 3:                    A3SS|ENSG00000158286|chrX|149608626-149608834||INC
#> 4:                    A3SS|ENSG00000158286|chrX|149608626-149608834||INC
#> 5: A3SS|ENSG00000158286|chrX|149608626-149608829|149608830-149608834|EXC
#> 6: A3SS|ENSG00000158286|chrX|149608626-149608829|149608830-149608834|EXC
#>    event_type event_id         gene_id    chr strand                 inc
#>        <char>   <char>          <char> <char> <char>              <char>
#> 1:       A3SS   A3SS:1 ENSG00000158286   chrX      - 149608626-149608834
#> 2:       A3SS   A3SS:1 ENSG00000158286   chrX      - 149608626-149608834
#> 3:       A3SS   A3SS:1 ENSG00000158286   chrX      - 149608626-149608834
#> 4:       A3SS   A3SS:1 ENSG00000158286   chrX      - 149608626-149608834
#> 5:       A3SS   A3SS:1 ENSG00000158286   chrX      - 149608626-149608829
#> 6:       A3SS   A3SS:1 ENSG00000158286   chrX      - 149608626-149608829
#>                    exc n_samples n_control n_case mean_psi_ctrl mean_psi_case
#>                 <char>     <int>     <int>  <int>         <num>         <num>
#> 1:                            NA        NA     NA            -1            -1
#> 2:                            NA        NA     NA            -1            -1
#> 3:                            NA        NA     NA            -1            -1
#> 4:                            NA        NA     NA            -1            -1
#> 5: 149608830-149608834        NA        NA     NA            -1            -1
#> 6: 149608830-149608834        NA        NA     NA            -1            -1
#>    delta_psi p.value  padj cooks_max   form     n n_used
#>        <num>   <num> <num>     <num> <char> <num>  <num>
#> 1:         1       0     0        -1    INC    -1     -1
#> 2:         1       0     0        -1    INC    -1     -1
#> 3:         1       0     0        -1    INC    -1     -1
#> 4:         1       0     0        -1    INC    -1     -1
#> 5:        -1       0     0        -1    EXC    -1     -1
#> 6:        -1       0     0        -1    EXC    -1     -1

Post-DI rMATS Import

Use get_rmats_post_di() to ingest already-computed rMATS result tables. rMATS reports IncLevelDifference as mean(IncLevel1) - mean(IncLevel2), where group 1 holds the --b1 samples. SpliceImpactR reads it as case - control, so by default rMATS group 1 is the case. If your case samples were --b2, set case_group = 2; with a file table you can also name the case group, for example case_group = "KO". Each import reports which group it treated as the case. In this example inclusion is higher in group 1, the case, so the inclusion form gets delta_psi = 1.

rmats_df <- data.frame(
  ID = 1L,
  GeneID = "ENSG00000182871",
  geneSymbol = "COL18A1",
  chr = "chr21",
  strand = "+",
  longExonStart_0base = 45505834L,
  longExonEnd = 45505966L,
  shortES = 45505837L,
  shortEE = 45505966L,
  flankingES = 45505357L,
  flankingEE = 45505431L,
  ID.2 = 2L,
  IJC_SAMPLE_1 = "5,5,5",
  SJC_SAMPLE_1 = "0,0,0",
  IJC_SAMPLE_2 = "0,0,0",
  SJC_SAMPLE_2 = "5,5,5",
  IncFormLen = 52L,
  SkipFormLen = 49L,
  PValue = 0.6967562,
  FDR = 1,
  IncLevel1 = "1.0,1.0,1.0",
  IncLevel2 = "0.0,0.0,0.0",
  IncLevelDifference = 1.0,
  stringsAsFactors = FALSE
)

user_res_rmats <- get_rmats_post_di(rmats_df, event_type = "A3SS")
#> get_rmats_post_di: rMATS group 1 (--b1) is the case; delta_psi = mean(IncLevel1) - mean(IncLevel2).
head(user_res_rmats)
#>    event_id event_type   form         gene_id    chr strand
#>      <char>     <char> <char>          <char> <char> <char>
#> 1:   A3SS:1       A3SS    INC ENSG00000182871  chr21      +
#> 2:   A3SS:1       A3SS    EXC ENSG00000182871  chr21      +
#>                                    inc               exc   p.value delta_psi
#>                                 <char>            <char>     <num>     <num>
#> 1: 45505358-45505431;45505835-45505966                   0.6967562         1
#> 2: 45505358-45505431;45505838-45505966 45505835-45505837 0.6967562        -1
#>     padj
#>    <num>
#> 1:     1
#> 2:     1

Transcript-Pair Entry Point

Use compare_transcript_pairs() when you want to start directly from chosen transcript pairs instead of event matching.

tx_ids <- unique(annotation_df$annotations$transcript_id)
tx_ids <- tx_ids[!is.na(tx_ids) & tx_ids != ""]

if (length(tx_ids) < 4L) {
  stop("Need at least four transcripts for transcript-pair example.")
}

transcript_pairs <- data.frame(
  transcript1 = tx_ids[1:2],
  transcript2 = tx_ids[3:4],
  stringsAsFactors = FALSE
)

user_matched <- compare_transcript_pairs(
  transcript_pairs = transcript_pairs,
  annotations = annotation_df$annotations
)
#> 2 out of 2 transcript pairs contained within annotations

head(user_matched)
#>    event_id event_type   form         gene_id    chr strand
#>      <char>     <char> <char>          <char> <char> <char>
#> 1:  TXCMP:1      TXCMP    EXC ENSG00000158286   chr1      +
#> 2:  TXCMP:1      TXCMP    INC ENSG00000158286   chr1      +
#> 3:  TXCMP:2      TXCMP    EXC ENSG00000158286   chr1      +
#> 4:  TXCMP:2      TXCMP    INC ENSG00000158286   chr1      +
#>                                                                inc
#>                                                             <char>
#> 1: 6205475-6206101;6206196-6206302;6206536-6206726;6207379-6208307
#> 2: 6205475-6206101;6206196-6206302;6206536-6206726;6207379-6208307
#> 3:                 6206129-6206302;6206536-6206726;6207379-6207564
#> 4:                 6206129-6206302;6206536-6206726;6207379-6207564
#>                                                                                                                                                                exc
#>                                                                                                                                                             <char>
#> 1:                                                                 6208618-6209025;6209115-6209196;6209268-6209539;6209924-6209970;6210223-6210295;6210370-6210900
#> 2:                                                                 6208618-6209025;6209115-6209196;6209268-6209539;6209924-6209970;6210223-6210295;6210370-6210900
#> 3: 6208860-6209025;6209414-6209539;6209924-6209970;6210223-6210295;6210370-6210438;6211867-6212053;6212218-6212416;6212682-6213183;6218289-6218369;6219236-6220805
#> 4: 6208860-6209025;6209414-6209539;6209924-6209970;6210223-6210295;6210370-6210438;6211867-6212053;6212218-6212416;6212682-6213183;6218289-6218369;6219236-6220805
#>    delta_psi p.value  padj n_samples n_control n_case   transcript_id
#>        <num>   <num> <num>     <num>     <num>  <num>          <char>
#> 1:        -1       0     0         0         0      0 ENST00000485539
#> 2:         1       0     0         0         0      0 ENST00000466994
#> 3:        -1       0     0         0         0      0 ENST00000496676
#> 4:         1       0     0         0         0      0 ENST00000484435
#>                                                                                                                                                              exons
#>                                                                                                                                                             <char>
#> 1:                                                                 ENSE00001934791;ENSE00003528881;ENSE00001948578;ENSE00003560800;ENSE00003516232;ENSE00001921263
#> 2:                                                                                                 ENSE00001843998;ENSE00001829482;ENSE00003686018;ENSE00001935169
#> 3: ENSE00001866172;ENSE00003662600;ENSE00003560800;ENSE00003516232;ENSE00003531305;ENSE00003488296;ENSE00001883823;ENSE00001858590;ENSE00003562079;ENSE00001844921
#> 4:                                                                                                                 ENSE00001551004;ENSE00003686018;ENSE00001828693

Session Info

sessionInfo()
#> R version 4.6.1 (2026-06-24)
#> Platform: x86_64-pc-linux-gnu
#> Running under: Ubuntu 26.04.1 LTS
#> 
#> Matrix products: default
#> BLAS:   /usr/lib/x86_64-linux-gnu/openblas-pthread/libblas.so.3 
#> LAPACK: /usr/lib/x86_64-linux-gnu/openblas-pthread/libopenblasp-r0.3.32.so;  LAPACK version 3.12.0
#> 
#> locale:
#>  [1] LC_CTYPE=en_US.UTF-8       LC_NUMERIC=C              
#>  [3] LC_TIME=en_US.UTF-8        LC_COLLATE=en_US.UTF-8    
#>  [5] LC_MONETARY=en_US.UTF-8    LC_MESSAGES=en_US.UTF-8   
#>  [7] LC_PAPER=en_US.UTF-8       LC_NAME=C                 
#>  [9] LC_ADDRESS=C               LC_TELEPHONE=C            
#> [11] LC_MEASUREMENT=en_US.UTF-8 LC_IDENTIFICATION=C       
#> 
#> time zone: Etc/UTC
#> tzcode source: system (glibc)
#> 
#> attached base packages:
#> [1] stats     graphics  grDevices utils     datasets  methods   base     
#> 
#> other attached packages:
#> [1] BiocParallel_1.47.0 SpliceImpactR_1.1.4 BiocStyle_2.41.0   
#> 
#> loaded via a namespace (and not attached):
#>   [1] splines_4.6.1               BiocIO_1.23.3              
#>   [3] bitops_1.1-0                ggplotify_0.1.3            
#>   [5] tibble_3.3.1                R.oo_1.27.1                
#>   [7] polyclip_1.10-7             enrichit_0.2.6             
#>   [9] XML_3.99-0.25               lifecycle_1.0.5            
#>  [11] httr2_1.3.0                 pwalign_1.9.1              
#>  [13] rstatix_1.1.0               processx_3.9.0             
#>  [15] lattice_0.23-1              MASS_7.3-66                
#>  [17] backports_1.5.1             magrittr_2.0.5             
#>  [19] sass_0.4.10                 rmarkdown_2.32             
#>  [21] jquerylib_0.1.4             yaml_2.3.12                
#>  [23] otel_0.2.0                  ggtangle_0.1.3             
#>  [25] cowplot_1.2.0               DBI_1.3.0                  
#>  [27] buildtools_1.0.0            RColorBrewer_1.1-3         
#>  [29] abind_1.4-8                 GenomicRanges_1.65.4       
#>  [31] purrr_1.2.2                 R.utils_2.13.0             
#>  [33] BiocGenerics_0.59.12        msigdbr_26.1.1             
#>  [35] RCurl_1.98-1.20             yulab.utils_0.2.5          
#>  [37] tweenr_2.0.3                rappdirs_0.3.4             
#>  [39] aisdk_1.4.12                gdtools_0.5.1              
#>  [41] IRanges_2.47.5              S4Vectors_0.51.10          
#>  [43] enrichplot_1.99.7           ggrepel_0.9.8              
#>  [45] tidytree_0.4.8              maketools_1.3.2            
#>  [47] DelayedArray_0.39.8         codetools_0.2-20           
#>  [49] DOSE_4.7.3                  ggforce_0.5.0              
#>  [51] tidyselect_1.2.1            aplot_0.3.2                
#>  [53] farver_2.1.2                matrixStats_1.5.0          
#>  [55] stats4_4.6.1                Seqinfo_1.3.2              
#>  [57] GenomicAlignments_1.49.2    jsonlite_2.0.0             
#>  [59] Formula_1.2-6               systemfonts_1.3.2          
#>  [61] tools_4.6.1                 ggnewscale_0.5.2           
#>  [63] treeio_1.37.1               PFAM.db_3.22.0             
#>  [65] Rcpp_1.1.2                  glue_1.8.1                 
#>  [67] SparseArray_1.13.4          gridExtra_2.3.1            
#>  [69] BiocBaseUtils_1.15.1        xfun_0.61                  
#>  [71] MatrixGenerics_1.25.0       qvalue_2.45.0              
#>  [73] dplyr_1.2.1                 withr_3.0.3                
#>  [75] BiocManager_1.30.27         fastmap_1.2.0              
#>  [77] callr_3.8.0                 digest_0.6.39              
#>  [79] R6_2.6.1                    gridGraphics_0.5-1         
#>  [81] GO.db_3.23.1                RSQLite_3.53.3             
#>  [83] cigarillo_1.3.1             R.methodsS3_1.8.2          
#>  [85] tidyr_1.3.2                 generics_0.1.4             
#>  [87] fontLiberation_0.1.0        data.table_1.18.6.1        
#>  [89] rtracklayer_1.73.0          S4Arrays_1.13.2            
#>  [91] httr_1.4.9                  htmlwidgets_1.6.4          
#>  [93] scatterpie_0.2.6            pkgconfig_2.0.3            
#>  [95] gtable_0.3.6                blob_1.3.0                 
#>  [97] S7_0.2.2                    XVector_0.53.0             
#>  [99] sys_3.4.3                   clusterProfiler_4.21.3     
#> [101] htmltools_0.5.9             fontBitstreamVera_0.1.1    
#> [103] carData_3.0-6               scales_1.4.0               
#> [105] Biobase_2.73.2              png_0.1-9                  
#> [107] ggfun_0.2.1                 knitr_1.52                 
#> [109] reshape2_1.4.5              rjson_0.2.23               
#> [111] nlme_3.1-171                curl_8.0.0                 
#> [113] org.Hs.eg.db_3.23.1         cachem_1.1.0               
#> [115] stringr_1.6.0               parallel_4.6.1             
#> [117] AnnotationDbi_1.75.2        restfulr_0.0.17            
#> [119] pillar_1.11.1               grid_4.6.1                 
#> [121] vctrs_0.7.3                 ggpubr_1.0.0               
#> [123] car_3.1-5                   tidydr_0.0.6               
#> [125] cluster_2.1.8.3             evaluate_1.0.5             
#> [127] cli_3.6.6                   compiler_4.6.1             
#> [129] Rsamtools_2.29.0            rlang_1.3.0                
#> [131] crayon_1.5.3                ggsignif_0.6.4             
#> [133] labeling_0.4.3              ps_1.9.3                   
#> [135] plyr_1.8.9                  fs_2.1.0                   
#> [137] ggiraph_0.9.6               stringi_1.8.9              
#> [139] assertthat_0.2.1            Biostrings_2.81.9          
#> [141] lazyeval_0.2.3              GOSemSim_2.39.3            
#> [143] fontquiver_0.2.1            Matrix_1.7-6               
#> [145] patchwork_1.3.2             bit64_4.8.6                
#> [147] ggplot2_4.0.3               KEGGREST_1.53.6            
#> [149] SummarizedExperiment_1.43.0 igraph_2.3.4               
#> [151] broom_1.0.13                memoise_2.0.1              
#> [153] bslib_0.12.0                ggtree_4.3.0               
#> [155] bit_4.6.0                   ape_5.8-1                  
#> [157] gson_0.2.1

SpliceImpactR