Introduction to IsoNetAnalyzeR

Introduction

A differential transcript usage (DTU) analysis reports isoforms whose share of their gene’s expression changes between conditions. When the isoforms of a gene encode different proteins, such a switch can change which interaction partners the gene product can reach. For each gene, IsoNetAnalyzeR takes the isoforms labelled up and the isoforms labelled down, expands each set into its neighbourhood in a predicted isoform-level interaction network (a virtual pulldown), and compares the two neighbourhoods.

The package provides:

  • per-gene comparison of up- and down-seeded neighbourhoods, with node, edge and topology statistics;
  • cross-gene summaries and plots;
  • connected components that join neighbourhoods across genes;
  • over-representation and gene set enrichment analysis;
  • diagnostics of the reference network;
  • resolution of versioned and protein-identical transcript identifiers.

It does not test for differential transcript usage, choose which isoforms are important, or measure interactions. Seed selection and its statistics belong to the upstream analysis. The network is a static prediction, so the comparisons describe predicted neighbourhoods, not measured gains or losses of interactions in either condition.

Installation

if (!requireNamespace("BiocManager", quietly = TRUE))
    install.packages("BiocManager")
BiocManager::install("IsoNetAnalyzeR")
library(IsoNetAnalyzeR)

The reference network

An analysis starts with a session object. Without further arguments it loads the bundled network: isoform-resolved predicted protein interactions produced by the Isoform Analysis Group at DTU. A Bayesian latent-variable model developed by Isabella Østerlund resolved high-confidence STRING physical interactions (combined score of at least 700; Szklarczyk et al. 2023) from genes to isoform pairs. Each vertex is a transcript; each edge carries the model’s interaction probability (pi) and is flagged as a predicted interaction (interact = 1). Provenance, including STRING’s CC BY 4.0 licence, is documented in system.file("PROVENANCE.md", package = "IsoNetAnalyzeR").

obj <- IsoNetAnalyzeR("GSE141432")
#> Loaded the bundled network: 58,164 transcripts, 158,458 edges.
#> Excluded 11,662 NMD and 69 NSD transcripts with their edges;
#> 46,433 transcripts and 146,267 edges remain.
obj
#> ── IsoNetAnalyzeR session ──────────────────────────────────────────────────────
#> Name        : GSE141432
#> Created     : 2026-11-03 19:50:07
#> Network     : 46,433 transcripts, 146,267 edges (bundled)
#> Analyses    : none
#> Lookups     : none
#> Diagnostics : not run
#> Metadata    : pi_threshold, network_cutoff, network_biotype_filter
#> ────────────────────────────────────────────────────────────────────────────────

New sessions exclude transcripts annotated (GENCODE v41) as nonsense-mediated decay (NMD) or non-stop decay (NSD) targets, both as seeds and as neighbours, because they are not expected to produce stable proteins. The exclusion is recorded in the session metadata:

eligibility <- get_metadata(obj, "network_biotype_filter")
unlist(eligibility[c("vertices_before", "excluded_nmd", "excluded_nsd",
                     "vertices_after", "edges_after")])
#> vertices_before    excluded_nmd    excluded_nsd  vertices_after     edges_after 
#>           58164           11662              69           46433          146267

Set include_nmd = TRUE or include_nsd = TRUE when creating the session to keep either class. A different network can be supplied as an igraph object through the network argument (Section @ref(isolated-and-unmapped-seeds) builds one).

The session also stores the neighbourhood cutoff calibrated when the network was built. During a pulldown, a neighbour is kept only if at least this fraction of its predicted interactions fall inside the candidate neighbourhood, which removes promiscuous hubs that would otherwise attach to every neighbourhood:

unname(get_metadata(obj, "network_cutoff"))
#> [1] 0.2608696

Input: isoform seeds

The input is a table with one row per seed isoform and three required columns:

Column Content
gene_id Gene identifier, matching the network’s gene attribute
isoform_id Transcript identifier; version suffixes are accepted
direction "up" or "down": labels for the two seed sets of a gene

The bundled network uses Ensembl gene and transcript identifiers.

The package ships seeds from GEO series GSE141432: white adipose tissue from 8 lean, 8 obese normoglycaemic and 8 obese type 2 diabetic (T2D) donors. The contrast is obese T2D versus lean. In the upstream analysis, genes passed a gene-level FDR of 0.05 and had a maximum absolute isoform fraction change (dIF) of at least 0.1. Within each gene, transcripts with an absolute dIF of at least 0.1 were tested with a Benjamini-Hochberg correction, stored here as within_gene_padj, and each gene kept at least one significant seed and seeds in both directions. Positive signed_dif values are labelled up.

seeds <- read.csv(system.file("extdata", "gse141432_t2d_vs_lean_seeds.csv",
                              package = "IsoNetAnalyzeR"))
head(seeds)
#>           gene_id gene_name      isoform_id direction signed_dif
#> 1 ENSG00000070371    CLTCL1 ENST00000621271        up   0.132414
#> 2 ENSG00000070371    CLTCL1 ENST00000427926      down  -0.151335
#> 3 ENSG00000182957   SPATA13 ENST00000424834      down  -0.202677
#> 4 ENSG00000182957   SPATA13 ENST00000382095        up   0.121520
#> 5 ENSG00000196074     SYCP2 ENST00000371001      down  -0.192732
#> 6 ENSG00000196074     SYCP2 ENST00000446834        up   0.191490
#>   within_gene_padj
#> 1      0.158062000
#> 2      0.027809500
#> 3      0.018969300
#> 4      0.049910500
#> 5      0.000715829
#> 6      0.000715829
table(seeds$direction)
#> 
#> down   up 
#>   70   67
length(unique(seeds$gene_id))
#> [1] 65

When preparing your own table, keep only the isoforms your study design selects, and label each by the sign of its usage change. A gene can have any number of seeds; a gene with seeds in only one direction is accepted but has nothing to be compared with.

Differential pulldown

run_differential_pulldown() expands every gene’s up and down seed sets by order interaction steps and applies the neighbourhood cutoff. With seed_padj_col, it also checks retained support: after identifier resolution, eligibility filtering and mapping, every gene must still have seeds in both directions and at least one seed with an adjusted p-value at or below seed_fdr. Complementary seeds do not have to be individually significant.

obj <- run_differential_pulldown(
  obj, "t2d", seeds,
  order = 1, confirm = FALSE, verbose = FALSE,
  seed_padj_col = "within_gene_padj", seed_fdr = 0.05
)
#> ── Network Coverage Report ─────────────────────────────────────────────────────
#> 
#> Isoforms
#>   Unique submitted IDs            : 137
#>   Unique resolved representatives : 137
#>   Duplicate-remapped IDs          : 0
#>   Version-normalized IDs          : 0
#>   Mapped, predicted edge          : 98
#>   Mapped, no predicted edge       : 39
#>   Unmapped                        : 0
#>   Mapping                         : 137 / 137 (100%)
#>   up: mapped with edge 58; mapped without edge 9; unmapped 0 (100% mapped)
#>   down: mapped with edge 40; mapped without edge 30; unmapped 0 (100% mapped)
#> 
#> Genes
#>   Total in input                  : 65
#>   At least 1 isoform in network   : 65 / 65 (100%)
#>   No submitted isoforms mapped    : 0
#>   Note: isoforms found in network with no interact==1 edges will produce a
#>         zero-interaction subgraph (single node), which is a valid result.
#> 
#> ────────────────────────────────────────────────────────────────────────────────
#> Done. Access the analysis with:
#>   get_analysis(obj, "t2d")

The coverage report distinguishes three outcomes for each seed: mapped with predicted interactions, mapped without any predicted interaction, and absent from the network. The tables are stored in the session:

coverage <- get_metadata(obj, "analysis_network_coverage")$t2d
coverage$by_direction[, c("direction", "n_input",
                          "n_mapped_with_predicted_edge",
                          "n_mapped_no_predicted_edge", "n_unmapped")]
#>   direction n_input n_mapped_with_predicted_edge n_mapped_no_predicted_edge
#> 1        up      67                           58                          9
#> 2      down      70                           40                         30
#>   n_unmapped
#> 1          0
#> 2          0
support <- get_metadata(obj, "analysis_seed_support")$t2d$gene_audit
table(support$reason)
#> 
#> retained 
#>       65

The cross-gene summary has one row per gene, with neighbourhood sizes and topology for each direction. Result tables are keyed by Ensembl gene ID; isonet_gene_names() returns the gene symbols of the bundled network for annotating them. Here are the genes with the largest difference in neighbourhood size:

gene_names <- isonet_gene_names()
cgs <- get_cross_gene_summary(obj, "t2d")
cgs$gene_name <- gene_names$gene_name[match(cgs$gene_id, gene_names$gene_id)]
cgs$delta_nodes <- cgs$up_nodes - cgs$down_nodes
head(cgs[order(-abs(cgs$delta_nodes)),
         c("gene_name", "up_nodes", "down_nodes", "up_edges", "down_edges",
           "jaccard")], 8)
#>                 gene_name up_nodes down_nodes up_edges down_edges jaccard
#> ENSG00000095002      MSH2       57          3      479          2  0.0000
#> ENSG00000091106     NLRC4       46          1      271          0  0.0000
#> ENSG00000015285       WAS       42          4      128          2  0.0222
#> ENSG00000142875    PRKACB        2         39        1        113  0.0250
#> ENSG00000141968      VAV1       40          4      153          3  0.0476
#> ENSG00000160460    SPTBN4        1         26        0        113  0.0000
#> ENSG00000080986     NDC80       23          1      111          0  0.0000
#> ENSG00000174996      KLC2       22          1       58          0  0.0000

jaccard is the overlap of the two neighbourhoods’ nodes: 1 means the isoforms reach the same predicted partners, 0 means entirely different ones.

plot_cross_gene_overview(obj, "t2d")
Neighbourhood statistics of every gene. Each panel plots the up-seeded pulldown (x) against the down-seeded one (y); genes on the diagonal have equal values.

Neighbourhood statistics of every gene. Each panel plots the up-seeded pulldown (x) against the down-seeded one (y); genes on the diagonal have equal values.

Example gene: CLTCL1

CLTCL1 encodes clathrin heavy chain-like 1 (CHC22). In this contrast, one CLTCL1 isoform gains usage and another loses it.

cltcl1 <- get_gene_result(obj, "t2d", "ENSG00000070371")
cltcl1
#> ── GeneDiffPulldown ────────────────────────────────────────────────────────────
#> Gene       : CLTCL1 (ENSG00000070371)
#> Up         : 1 isoform(s) input → subgraph with 8 nodes, 8 edges
#> Down       : 1 isoform(s) input → subgraph with 12 nodes, 33 edges
#> ────────────────────────────────────────────────────────────────────────────────
stats <- get_isoform_stats(obj, "t2d", "ENSG00000070371")
stats[stats$is_seed, c("direction", "isoform", "degree", "mean_pi")]
#>    direction         isoform degree mean_pi
#> 8         up ENST00000621271      7  0.2005
#> 20      down ENST00000427926     11  0.3185

The up-labelled isoform reaches 7 predicted partners, transcripts of AMPH, AP1B1, BIN1, EPN2 and EPN3, through 8 edges. The down-labelled isoform’s neighbourhood has 11 partners, transcripts of AP1B1, AP2A1, AP2A2, CLTA, CLTB, DNAJC6, EHD2, HIP1 and PICALM, and 33 edges. The overlay shows both neighbourhoods at once, with each gene’s transcripts merged into one node:

plot_network(cltcl1, direction = "overlay", level = "gene",
             node_size = 3.5, label_size = 3)
Predicted neighbourhoods of the up- and down-labelled CLTCL1 isoforms, one node per gene. The CLTCL1 node (diamond) holds both seeds. Genes and edges are green when reached only from the up-labelled isoform, orange with dashed edges when reached only from the down-labelled isoform, and blue when reached from both.

Predicted neighbourhoods of the up- and down-labelled CLTCL1 isoforms, one node per gene. The CLTCL1 node (diamond) holds both seeds. Genes and edges are green when reached only from the up-labelled isoform, orange with dashed edges when reached only from the down-labelled isoform, and blue when reached from both.

Nodes are labelled with gene symbols. The default, level = "isoform", draws one node per transcript instead, so a gene with several isoforms in a neighbourhood appears more than once and the two seeds are separate triangles. To label nodes with transcript identifiers, or to hide the labels, use the node_labels argument, which plot_component() also accepts:

plot_network(cltcl1, direction = "overlay")                       # isoforms
plot_network(cltcl1, direction = "overlay", node_labels = "id")   # IDs
plot_network(cltcl1, direction = "overlay", node_labels = "none") # no labels

Read this as a statement about the reference network: the isoform used more in obese T2D tissue is predicted to reach fewer partners than the isoform used less, and different ones. Whether those partners are bound in adipose tissue needs independent evidence.

Connected components

Neighbourhoods of different genes can share nodes. Component analysis merges all neighbourhoods into one graph and reports its connected components, so that genes whose isoform switches converge on the same proteins appear together.

obj <- run_component_analysis(obj, "t2d", label = "comp")
#> Building merged cross-gene graph for 't2d' (direction = 'both', level = 'isoform') ...
#> Decomposing into connected components ...
#> 122 component(s) retained.
#> Computing Tier-1 structural statistics ...
#> Computing Tier-2 topological statistics ...
#>   672 node(s) total  |  10 mixed component(s)  |  13 shared hub(s) scored.
#> Done. Access the results with:
#>   get_component_overview(obj, "t2d", "comp")
overview <- get_component_overview(obj, "t2d", "comp")
head(overview[, c("component_id", "n_nodes", "n_edges", "n_seeds_up",
                  "n_seeds_down", "mixed", "n_shared_nodes",
                  "topology_class")])
#>   component_id n_nodes n_edges n_seeds_up n_seeds_down mixed n_shared_nodes
#> 1            1      89     308          2            3  TRUE              6
#> 2            2      57     479          1            0 FALSE              0
#> 3            3      46     271          1            0 FALSE              0
#> 4            4      40     114          1            1  TRUE              0
#> 5            5      33     136          0            2 FALSE              7
#> 6            6      23     111          1            0 FALSE              0
#>   topology_class
#> 1     chain-like
#> 2      hub-spoke
#> 3      hub-spoke
#> 4      hub-spoke
#> 5     chain-like
#> 6     chain-like

A mixed component contains both up and down seeds. topology_class is a heuristic description of shape (hub-spoke, chain-like or clique-like), not a mechanistic classification.

The largest mixed component joins the neighbourhoods of GAB2, VAV1 and WAS. With colour_neighbours = TRUE, each neighbour is coloured by the seed direction it is closest to:

plot_component(obj, "t2d", "comp", component_id = largest_mixed,
               colour_by = "direction", colour_neighbours = TRUE)
The largest mixed component. Diamonds are seeds. Every other node is coloured by the seed direction it is closest to; purple nodes are equally close to up and down seeds.

The largest mixed component. Diamonds are seeds. Every other node is coloured by the seed direction it is closest to; purple nodes are equally close to up and down seeds.

Functional annotation

Over-representation in components

run_component_ora() tests each component’s genes against a gene set collection, with all genes in the network as the universe. Without a pathways argument it retrieves GO gene sets through msigdbr; for a reproducible analysis, pass a saved collection. Here we use the MSigDB Hallmark collection shipped with the package (50 gene sets, Ensembl IDs):

hallmark <- readRDS(system.file("extdata", "msigdb_hallmark_ensembl.rds",
                                package = "IsoNetAnalyzeR"))
attr(hallmark, "source")$db_version
#> [1] "2026.1.Hs"
obj <- run_component_ora(obj, "t2d", "comp", label = "hallmark",
                         pathways = hallmark)
#> Building isoform-to-gene mapping from network vertex attributes ...
#> Universe: 10530 unique gene(s) represented in the full network.
#> Using user-supplied pathway list (50 gene set(s)).
#> Extracting gene sets per component (nodes = 'all') ...
#> Running fgsea::fora() on 122 component(s) ...
#> Done. 122 component(s) tested, 0 skipped. 30 significant term(s) at padj < 0.05.
#> Access results with:
#>   get_component_ora_result(obj, "t2d", "comp", "hallmark")
ora <- get_ora_table(obj, "t2d", "comp", "hallmark", padj_cutoff = 0.05,
                     top_n = 1)
head(ora[, c("component_id", "pathway", "overlap", "padj")], 8)
#>   component_id                            pathway overlap     padj
#> 1            1           HALLMARK_APICAL_JUNCTION      10 0.000001
#> 2            2                HALLMARK_DNA_REPAIR       7 0.000006
#> 3            3 HALLMARK_INTERFERON_GAMMA_RESPONSE       5 0.000853
#> 4            5           HALLMARK_MITOTIC_SPINDLE       3 0.041556
#> 5            6            HALLMARK_G2M_CHECKPOINT       5 0.000716
#> 6            7       HALLMARK_ALLOGRAFT_REJECTION       8 0.000000
#> 7            9           HALLMARK_MITOTIC_SPINDLE       5 0.000038
#> 8           11       HALLMARK_ALLOGRAFT_REJECTION       5 0.000058

Gene set enrichment across genes

run_fgsea() ranks the analysed genes by a neighbourhood statistic, by default the difference in node count between the up and down neighbourhoods, and runs fgsea. With 65 genes the ranked list is short, so only small gene sets can be tested and power is low. We lower minSize for illustration.

set.seed(42)
obj <- run_fgsea(obj, "t2d", label = "hallmark_gsea", pathways = hallmark,
                 minSize = 3)
#> Ranking genes by preset 'delta_nodes'  [up_nodes - down_nodes].
#> Using user-supplied pathway list (50 gene sets).
#> Running fgsea: 65 genes ranked, 50 pathway(s), minSize = 3, maxSize = 500, nperm = 1000 ...
#> Done. 3 pathway(s) tested; 0 significant at padj < 0.05.
#> Access results with:
#>   get_gsea_result(obj, "t2d", "hallmark_gsea")
get_gsea_table(obj, "t2d", "hallmark_gsea")[, c("pathway", "size", "NES",
                                                "pval", "padj")]
#>                        pathway size        NES      pval      padj
#> 1 HALLMARK_ALLOGRAFT_REJECTION    3  1.1775949 0.2156573 0.4656652
#> 2   HALLMARK_KRAS_SIGNALING_UP    4  1.1292622 0.3104435 0.4656652
#> 3     HALLMARK_HEME_METABOLISM    3 -0.6924581 0.8098160 0.8098160

Custom rankings can be passed through the stats argument.

Network diagnostics

Results depend on the structure of the reference network. Diagnostics report its components, degree distribution, clustering and sticky nodes, high-degree vertices that join many neighbourhoods. This step takes about half a minute on the bundled network.

obj <- run_network_diagnostics(obj)
#> Running network diagnostics...
#>   Building predicted subgraphs...
#>   Computing connected components...
#>   Computing degree distributions...
#>   Computing clustering coefficients...
#>   Computing per-gene connectivity status...
#>   Computing global topology metrics...
#>   Identifying sticky nodes...
#> Network diagnostics complete. Access with get_network_diagnostics(obj).
get_network_diagnostics(obj)
#> ── NetworkDiagnostics ──────────────────────────────────────────────────────────
#> Created    : 2026-11-03 19:51:51
#> 
#> Full Network
#>   Isoforms / genes     : 46433 / 10530
#>   Total edges          : 146267
#>   Predicted edges      : 146267
#>   Isolated isoforms    : 25617 (55.2%)
#> 
#> Predicted Isoform Network
#>   Isoforms / genes     : 20816 / 10458
#>   Edges                : 146267
#>   Density              : 7e-04
#>   Avg degree           : 14.05
#>   Avg clustering       : 0.4273
#>   Components           : 438
#>   Largest component    : 19352 nodes (93%)
#>   Diameter             : 17
#>   Avg path length      : 5.1
#> 
#> Predicted Gene-Level Network
#>   Genes / edges        : 10458 / 81004
#> 
#> Sticky Nodes  [quantile / cutoff = 0.95]
#>   Sticky isoforms      : 1019 (threshold degree > 58)
#>   Sticky genes         : 511 (threshold degree > 70)
#> 
#> Gene Connectivity Status
#>   All connected        : 3070 genes
#>   Partial              : 7388 genes
#>   All isolated         : 72 genes
#> ────────────────────────────────────────────────────────────────────────────────
plot_network_diagnostics(obj)
Diagnostics of the bundled network after decay-transcript exclusion.

Diagnostics of the bundled network after decay-transcript exclusion.

Single-isoform lookups

To explore one isoform without a differential design, use lookup_virtual_pulldown(). Seeds can be isoforms or, with level = "gene", all of a gene’s isoforms in the network.

obj <- lookup_virtual_pulldown(obj, "cltcl1_down", seed = "ENST00000427926",
                               order = 1)
#> Using stored network_cutoff: 0.2609
#> Running virtual pulldown: 1 seed vertex/vertices, order = 1, cutoff = 0.2609
#> Look-up 'cltcl1_down' stored. Access with: get_lookup(obj, "cltcl1_down")
lookup <- get_lookup_stats(obj, "cltcl1_down")
subgraph <- get_lookup_subgraph(obj, "cltcl1_down")
gene_of <- stats::setNames(igraph::V(subgraph)$gene, igraph::V(subgraph)$name)
lookup$gene_name <- gene_names$gene_name[
  match(gene_of[lookup$isoform], gene_names$gene_id)
]
rownames(lookup) <- NULL
head(lookup[, c("isoform", "gene_name", "is_seed", "degree",
                "distance_from_seed")])
#>           isoform gene_name is_seed degree distance_from_seed
#> 1 ENST00000345519      CLTA   FALSE      8                  1
#> 2 ENST00000263277      EHD2   FALSE      5                  1
#> 3 ENST00000345807      CLTB   FALSE      8                  1
#> 4 ENST00000448903     AP2A2   FALSE      7                  1
#> 5 ENST00000336926      HIP1   FALSE      4                  1
#> 6 ENST00000357586     AP1B1   FALSE      3                  1

Identical proteins

Different transcripts can encode the same protein. The package ships a map, supplied by Isabella Østerlund, of 16,630 redundant transcripts to their protein-identical representatives. Seed identifiers are resolved against it by default: version suffixes are removed, and a redundant transcript is replaced by its representative when the network and gene annotation agree. If a gene’s up and down seeds resolve to the same representative, the apparent switch is between identical proteins, so that gene is dropped.

map <- isonet_duplicate_map()
ids <- c("ENST00000427926.5", map$redundant_ENST[map$network_applicable][1])
resolved <- resolve_duplicate_isoforms(ids, network = get_network(obj))
attr(resolved, "dedup_audit")[, c("submitted_id", "resolved_id",
                                  "resolution")]
#>        submitted_id     resolved_id         resolution
#> 1 ENST00000427926.5 ENST00000427926 version_normalized
#> 2   ENST00000360466 ENST00000349431 duplicate_remapped

For DTU testing itself, aggregate counts of protein-identical transcripts before fitting models, so that isoform fractions are computed over distinct proteins. Gene and library totals are preserved:

counts <- matrix(
  c(10L, 4L, 7L,
    2L, 1L, 3L,
    5L, 5L, 5L),
  nrow = 3, byrow = TRUE,
  dimnames = list(c(map$kept_ENST[1], map$redundant_ENST[1], map$kept_ENST[2]),
                  paste0("sample", seq_len(3)))
)
features <- data.frame(
  gene_id = c(map$kept_ENSG[1], map$redundant_ENSG[1], map$kept_ENSG[2]),
  isoform_id = rownames(counts)
)
aggregated <- aggregate_duplicate_counts(counts, features)
aggregated$object
#>                 sample1 sample2 sample3
#> ENST00000433179      12       5      10
#> ENST00000421241       5       5       5

Isolated vs unmapped

A seed in the network without any predicted interaction gives a valid one-node neighbourhood. A seed absent from the network gives no neighbourhood at all. The two must not be confused, and neither shows that the isoform lacks interactions. This small custom network has one transcript with neighbours, one without, and a seed that is missing:

vertices <- data.frame(
  name = c("tx_A_up", "tx_A_down", "tx_B", "tx_C"),
  gene = c("gene_A", "gene_A", "gene_B", "gene_C"),
  transcript_biotype = "protein_coding"
)
edges <- data.frame(from = c("tx_A_up", "tx_A_down"), to = "tx_B",
                    interact = 1L, pi = c(0.90, 0.75))
toy <- IsoNetAnalyzeR(
  "toy",
  network = igraph::graph_from_data_frame(edges, directed = FALSE,
                                          vertices = vertices),
  metadata = list(network_cutoff = 0)
)
#> Using the supplied network: 4 transcripts, 2 edges.
#> Excluded 0 NMD and 0 NSD transcripts with their edges;
#> 4 transcripts and 2 edges remain.
toy_seeds <- data.frame(
  gene_id = c("gene_A", "gene_A", "gene_C", "gene_D"),
  isoform_id = c("tx_A_up", "tx_A_down", "tx_C", "tx_absent"),
  direction = c("up", "down", "up", "down")
)
toy <- run_differential_pulldown(toy, "toy", toy_seeds, order = 1,
                                 confirm = FALSE, verbose = FALSE)
#> ── Network Coverage Report ─────────────────────────────────────────────────────
#> 
#> Isoforms
#>   Unique submitted IDs            : 4
#>   Unique resolved representatives : 4
#>   Duplicate-remapped IDs          : 0
#>   Version-normalized IDs          : 0
#>   Mapped, predicted edge          : 2
#>   Mapped, no predicted edge       : 1
#>   Unmapped                        : 1
#>   Mapping                         : 3 / 4 (75%)
#>   up: mapped with edge 1; mapped without edge 1; unmapped 0 (100% mapped)
#>   down: mapped with edge 1; mapped without edge 0; unmapped 1 (50% mapped)
#> 
#> Genes
#>   Total in input                  : 3
#>   At least 1 isoform in network   : 2 / 3 (66.7%)
#>   No submitted isoforms mapped    : 1
#>   Note: isoforms found in network with no interact==1 edges will produce a
#>         zero-interaction subgraph (single node), which is a valid result.
#> 
#> ────────────────────────────────────────────────────────────────────────────────
#> Warning: Gene 'gene_D' (down): virtual_pulldown failed - None of the seed
#> vertices were found in the network.
#> Done. Access the analysis with:
#>   get_analysis(toy, "toy")
isolated <- get_pulldown_graph(get_gene_result(toy, "toy", "gene_C"), "up")
c(nodes = igraph::vcount(isolated), edges = igraph::ecount(isolated))
#> nodes edges 
#>     1     0
is.null(get_pulldown_graph(get_gene_result(toy, "toy", "gene_D"), "down"))
#> [1] TRUE

The run warns that gene_D has no seed in the network.

Interpreting results

Both neighbourhoods of a gene come from the same static network, so:

  • node and edge differences compare seed-derived neighbourhoods, not condition-specific interactomes;
  • n_edges_gained and n_edges_lost in the cross-gene summary count edges unique to the up-seeded and down-seeded neighbourhoods respectively;
  • a mixed component only records that both seed labels are present;
  • a high-degree or high-betweenness node is a structural observation that needs biological validation.

Acknowledgements

Isabella Østerlund developed the Bayesian model behind the isoform-resolved interaction network and supplied the map of protein-identical transcripts. She and Kristoffer Vitting-Seerup (Isoform Analysis Group, DTU) supervised the thesis project in which the package was first developed. The network derives from STRING (Szklarczyk et al. 2023), used under the CC BY 4.0 licence.

Session information

sessionInfo()
#> R version 4.6.2 (2026-10-29)
#> 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] IsoNetAnalyzeR_0.99.0 BiocStyle_2.41.0     
#> 
#> loaded via a namespace (and not attached):
#>  [1] fastmatch_1.1-8     gtable_0.3.6        xfun_0.62          
#>  [4] bslib_0.12.0        ggplot2_4.0.3       ggrepel_0.9.8      
#>  [7] lattice_0.23-1      vctrs_0.7.3         tools_4.6.2        
#> [10] generics_0.1.4      parallel_4.6.2      tibble_3.3.1       
#> [13] pkgconfig_2.0.3     Matrix_1.7-6        data.table_1.18.6.1
#> [16] RColorBrewer_1.1-3  S7_0.2.2            lifecycle_1.0.5    
#> [19] compiler_4.6.2      farver_2.1.2        ggforce_0.5.0      
#> [22] fgsea_1.39.2        graphlayouts_1.2.5  codetools_0.2-20   
#> [25] htmltools_0.5.9     sys_3.4.3           buildtools_1.0.0   
#> [28] sass_0.4.10         yaml_2.3.12         pillar_1.11.1      
#> [31] jquerylib_0.1.4     tidyr_1.3.2         MASS_7.3-66        
#> [34] BiocParallel_1.47.0 cachem_1.1.0        viridis_0.6.5      
#> [37] nlme_3.1-171        tidyselect_1.2.1    digest_0.6.39      
#> [40] dplyr_1.2.1         purrr_1.2.2         labeling_0.4.3     
#> [43] maketools_1.3.2     splines_4.6.2       cowplot_1.2.0      
#> [46] polyclip_1.10-7     fastmap_1.2.0       grid_4.6.2         
#> [49] cli_3.6.6           magrittr_2.0.5      patchwork_1.3.2    
#> [52] ggraph_2.2.2        tidygraph_1.3.1     withr_3.0.3        
#> [55] scales_1.4.0        rmarkdown_2.32      igraph_2.3.4       
#> [58] otel_0.2.0          gridExtra_2.3.1     progressr_1.1.0    
#> [61] memoise_2.0.1       evaluate_1.0.5      knitr_1.52         
#> [64] viridisLite_0.4.3   mgcv_1.9-4          rlang_1.3.0        
#> [67] Rcpp_1.1.2          glue_1.8.1          tweenr_2.0.3       
#> [70] BiocManager_1.30.27 jsonlite_2.0.0      R6_2.6.1