SpliceImpactR
October 9, 2026
Package
SpliceImpactR 1.1.4
Contents
- 1 Overview
- 2 Load Reference Resources
- 3 Build Sample Manifest
- 4 Read Splicing Events
- 5 Sample-Level QC and Exploration
- 6 Differential Inclusion (DI)
- 7 Match Events to Sequences and Build Transcript Pairs
- 8 Sequence-Level Consequences
- 9 Domain-Level Changes and Domain Enrichment
- 10 PPI Rewiring
- 11 Gene Enrichment Foreground Probes
- 12 Visualize Transcript and Protein Consequences
- 13 Inspect a Single Event
- 14 Integrative Visualization
- 15 Output Description (
hits_final,data,res) - 16 S4 Applications and Accessors
- 17 Custom Input Entry Points
- 18 Session Info
1 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.
1.1 Installation
Install with one of the following methods.
1.1.1 Option 1: BiocManager
if (!requireNamespace("BiocManager", quietly = TRUE)) {
install.packages("BiocManager")
}
BiocManager::install("SpliceImpactR")1.1.2 Option 2: devtools (GitHub)
if (!requireNamespace("devtools", quietly = TRUE)) {
install.packages("devtools")
}
devtools::install_github("fiszbein-lab/SpliceImpactR")1.2 External tools and acronyms
rMATS(replicate Multivariate Analysis of Transcript Splicing) provides event-level splicing tables.HITindexquantifies first/last exon usage and labels hybrid exons.PPIDMProtein-Protein Interaction Domain Miner for domain-domain Interaction derived from PPI and 3did’s domain-domain interactionsELMEukaryotic Linear Motif Database for accessing Short Linear Motif occurences and domain-motif interactionsBiomaRtaccesses data from InterPro, PFAM, SignalP, TMHMM, CDD, Mobidb-lite
1.3 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.
2 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 133 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
)3.1 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.
4 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 41If 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 415 Sample-Level QC and Exploration
5.1 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
Figure 1: 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 15.2 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
Figure 2: Depth-normalized event counts and PSI overview
6 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 72plot_di_volcano_dt(res)
Figure 3: Differential inclusion volcano plot
7 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_psirows are the case and negative rows the control, for eachevent_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 thanmax_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 = FALSEdisables this). Comparisons still without a pair are listed inunmatchedwith 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 16With 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.
7.1 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, andgene_idto 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 positivedelta_psirows (case) with all negative rows (control) for eachevent_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
Figure 4: Proximal versus distal terminal exon summary
8 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
Figure 5: Protein alignment summary by event type
length_plot <- plot_length_comparison(seq_compare)
length_plot
Figure 6: Case versus control isoform length comparison
9 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 10enriched_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
Figure 7: 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] 110 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
Figure 8: Summary of predicted case/control PPI changes
11 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_finalorhits_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.
12 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
Figure 9: Transcript and protein views for one paired event
protein_centric
Figure 10: Transcript and protein views for one paired event
13 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
Figure 11: 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 INC14 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
Figure 12: 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.0000000015 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.
15.1 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_absderivatives. - 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), andsummary_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.
15.2 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.
15.3 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.
16 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 33obj_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 1fg_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"16.1 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)17 Custom Input Entry Points
These helpers support non-standard inputs at multiple stages of the workflow.
17.1 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-836184117.2 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: case17.3 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 -117.4 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: 117.5 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;ENSE0000182869318 Session Info
sessionInfo()
#> R version 4.6.1 (2026-06-24)
#> Platform: x86_64-pc-linux-gnu
#> Running under: Ubuntu 24.04.5 LTS
#>
#> Matrix products: default
#> BLAS: /home/biocbuild/bbs-3.24-bioc/R/lib/libRblas.so
#> LAPACK: /usr/lib/x86_64-linux-gnu/lapack/liblapack.so.3.12.0 LAPACK version 3.12.0
#>
#> locale:
#> [1] LC_CTYPE=en_US.UTF-8 LC_NUMERIC=C
#> [3] LC_TIME=en_GB LC_COLLATE=C
#> [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: America/New_York
#> 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] RColorBrewer_1.1-3 abind_1.4-8
#> [29] GenomicRanges_1.65.4 purrr_1.2.2
#> [31] R.utils_2.13.0 BiocGenerics_0.59.12
#> [33] msigdbr_26.1.1 RCurl_1.98-1.20
#> [35] yulab.utils_0.2.5 tweenr_2.0.3
#> [37] rappdirs_0.3.4 aisdk_1.4.12
#> [39] gdtools_0.5.1 IRanges_2.47.5
#> [41] S4Vectors_0.51.10 enrichplot_1.99.8
#> [43] ggrepel_0.9.8 tidytree_0.4.8
#> [45] DelayedArray_0.39.8 codetools_0.2-20
#> [47] DOSE_4.7.3 ggforce_0.5.0
#> [49] tidyselect_1.2.1 aplot_0.3.2
#> [51] farver_2.1.2 matrixStats_1.5.0
#> [53] stats4_4.6.1 Seqinfo_1.3.2
#> [55] GenomicAlignments_1.49.2 jsonlite_2.0.0
#> [57] Formula_1.2-6 systemfonts_1.3.2
#> [59] tools_4.6.1 ggnewscale_0.5.2
#> [61] treeio_1.37.1 PFAM.db_3.22.0
#> [63] Rcpp_1.1.2 glue_1.8.1
#> [65] SparseArray_1.13.4 gridExtra_2.3.1
#> [67] BiocBaseUtils_1.15.1 xfun_0.61
#> [69] MatrixGenerics_1.25.0 qvalue_2.45.0
#> [71] dplyr_1.2.1 withr_3.0.3
#> [73] BiocManager_1.30.27 fastmap_1.2.0
#> [75] callr_3.8.0 digest_0.6.39
#> [77] R6_2.6.1 gridGraphics_0.5-1
#> [79] GO.db_3.23.1 dichromat_2.0-1
#> [81] RSQLite_3.53.3 cigarillo_1.3.1
#> [83] R.methodsS3_1.8.2 tidyr_1.3.2
#> [85] generics_0.1.4 fontLiberation_0.1.0
#> [87] data.table_1.18.6.1 rtracklayer_1.73.0
#> [89] S4Arrays_1.13.2 httr_1.4.9
#> [91] htmlwidgets_1.6.4 scatterpie_0.2.6
#> [93] pkgconfig_2.0.3 gtable_0.3.6
#> [95] blob_1.3.0 S7_0.2.2
#> [97] XVector_0.53.0 clusterProfiler_4.21.3
#> [99] htmltools_0.5.9 fontBitstreamVera_0.1.1
#> [101] carData_3.0-6 bookdown_0.48
#> [103] scales_1.4.0 Biobase_2.73.2
#> [105] png_0.1-9 ggfun_0.2.2
#> [107] knitr_1.52 reshape2_1.4.5
#> [109] rjson_0.2.23 nlme_3.1-171
#> [111] curl_8.0.0 org.Hs.eg.db_3.23.1
#> [113] cachem_1.1.0 stringr_1.6.0
#> [115] parallel_4.6.1 AnnotationDbi_1.75.2
#> [117] restfulr_0.0.17 pillar_1.11.1
#> [119] grid_4.6.1 vctrs_0.7.3
#> [121] ggpubr_1.0.0 car_3.1-5
#> [123] tidydr_0.0.7 cluster_2.1.8.3
#> [125] evaluate_1.0.5 tinytex_0.61
#> [127] magick_2.9.1 cli_3.6.6
#> [129] compiler_4.6.1 Rsamtools_2.29.0
#> [131] rlang_1.3.0 crayon_1.5.3
#> [133] ggsignif_0.6.4 labeling_0.4.3
#> [135] ps_1.9.3 plyr_1.8.9
#> [137] fs_2.1.0 ggiraph_0.9.6
#> [139] stringi_1.8.9 assertthat_0.2.1
#> [141] Biostrings_2.81.9 lazyeval_0.2.3
#> [143] GOSemSim_2.39.3 fontquiver_0.2.1
#> [145] Matrix_1.7-6 patchwork_1.3.2
#> [147] bit64_4.8.6 ggplot2_4.0.3
#> [149] KEGGREST_1.53.6 SummarizedExperiment_1.43.0
#> [151] igraph_2.3.4 broom_1.0.13
#> [153] memoise_2.0.1 bslib_0.12.0
#> [155] ggtree_4.3.0 bit_4.6.0
#> [157] ape_5.8-1 gson_0.2.2