Conserved and Cell Line-Specific Signaling in Sézary Syndrome

1 Introduction

Cell–cell communication predictions from LIANA (L1–L7) were validated against curated ligand–receptor pairs from the OmniPath LigRecExtra database (downloaded 2026-03-13 from https://omnipathdb.org/interactions?datasets=ligrecextra).

---

2 Setup

---

3 Load LIANA Results

read_liana <- function(path, name) {
  df <- read.csv(path)
  df$CellLine <- name
  return(df)
}

L1 <- read_liana("liana_L1_aggregate_results.csv", "L1")
L2 <- read_liana("liana_L2_aggregate_results.csv", "L2")
L3 <- read_liana("liana_L3_aggregate_results.csv", "L3")
L4 <- read_liana("liana_L4_aggregate_results.csv", "L4")
L5 <- read_liana("liana_L5_aggregate_results.csv", "L5")
L6 <- read_liana("liana_L6_aggregate_results.csv", "L6")
L7 <- read_liana("liana_L7_aggregate_results.csv", "L7")

all_data <- bind_rows(L1, L2, L3, L4, L5, L6, L7)

# Create unique interaction label
all_data$Interaction <- paste(
  all_data$ligand.complex, "->", all_data$receptor.complex
)

cat("Total LIANA records loaded:", nrow(all_data), "\n")

Total LIANA records loaded: 48958 

cat("Unique LR labels:          ", length(unique(all_data$Interaction)), "\n")

Unique LR labels:           605 

---

4 Load OmniPath LR Reference (Static CSV — No Package Needed)

# SOLUTION: Load from pre-downloaded static file
# This completely bypasses the OmnipathR package and the httr2 incompatibility
# The file was downloaded once from:
# https://omnipathdb.org/interactions?datasets=ligrecextra&organisms=9606&genesymbols=yes

lr_db <- read.csv("OmniPath_LigRecExtra_Human_2026.csv",
                  stringsAsFactors = FALSE)

# Validate it loaded correctly
stopifnot(all(c("ligand", "receptor") %in% colnames(lr_db)))

# Create a fast lookup key for matching
lr_db <- lr_db %>%
  dplyr::distinct(ligand, receptor) %>%
  dplyr::filter(nchar(ligand) > 0 & nchar(receptor) > 0) %>%
  dplyr::mutate(lr_key = paste(ligand, receptor, sep = "__"))

cat("OmniPath LR reference pairs loaded:", nrow(lr_db), "\n")

OmniPath LR reference pairs loaded: 6555 

---

5 Filter LIANA Predictions Against OmniPath

# Create matching key in LIANA data
all_data <- all_data %>%
  dplyr::mutate(lr_key = paste(ligand.complex, receptor.complex, sep = "__"))

# Semi-join: retain only LIANA pairs confirmed in the OmniPath LR database
filtered_data <- all_data %>%
  dplyr::filter(lr_key %in% lr_db$lr_key)

cat("Before OmniPath filter:", nrow(all_data),     "\n")

Before OmniPath filter: 48958 

cat("After OmniPath filter: ", nrow(filtered_data), "\n")

After OmniPath filter:  31427 

cat("Retention rate:        ",
    round(nrow(filtered_data) / nrow(all_data) * 100, 1), "%\n")

Retention rate:         64.2 %

---

6 Remove Biologically Implausible Interactions

Manually exclude known false positives: intracellular adaptor proteins listed as ligands/receptors, and Sézary-specific invalid pairs (e.g., DPP4/CD26, which is absent on Sézary cells by definition).

# Intracellular proteins that cannot mediate extracellular signaling
intracellular_blacklist <- c(
  "FADD",  # Death domain adaptor — intracellular
  "TRADD", # Death domain adaptor — intracellular
  "TRAF2", # Cytoplasmic TNFR adaptor
  "GNAI2", # G-protein alpha subunit — cytoplasmic
  "PTGS2"  # COX-2 enzyme — intracellular
)

# Sézary-specific invalid receptors (CD26/DPP4 is absent on malignant Sézary cells)
sezary_blacklist_receptors <- c("DPP4")

filtered_data <- filtered_data %>%
  dplyr::filter(
    !ligand.complex   %in% intracellular_blacklist,
    !receptor.complex %in% intracellular_blacklist,
    !receptor.complex %in% sezary_blacklist_receptors
  )

cat("After biological plausibility filter:", nrow(filtered_data), "\n")

After biological plausibility filter: 29399 

# ─── FINAL SAFETY-VERIFIED BIOLOGICAL FILTER ─────────────────────────────

# CATEGORY 1: Confirmed intracellular — safe to blacklist
# All verified: expressed in T cells but NO extracellular domain or secreted form
intracellular_blacklist <- c(
  "FADD",    # Cytoplasmic death domain adaptor [web:44]
  "TRADD",   # Cytoplasmic TNFR1 adaptor [web:44]
  "TRAF2",   # Cytoplasmic TNFR adaptor [web:106][web:115]
  "TRAF5",   # Cytoplasmic IL-6R adaptor in CD4 T cells [web:109] — intracellular only
  "TRAF6",   # Cytoplasmic TCR/CD28 adaptor [web:106]
  "GNAI2",   # G-protein alpha — no secreted isoform
  "GNAI3",   # G-protein alpha — no secreted isoform
  "PTGS2",   # COX-2 — ER membrane enzyme, not secreted
  "RIPK1",   # Cytoplasmic kinase [web:110][web:113]
  "CASP8",   # Cytoplasmic caspase [web:110]
  "BIRC3",   # Cytoplasmic IAP
  "MAP3K7"   # TAK1 kinase — cytoplasmic [web:119]
)

# CATEGORY 2: Neuronal — no evidence of expression in any T cell subset
neuronal_blacklist <- c(
  "NLGN1", "NLGN2", "NLGN3",   # Neuroligins [web:71][web:77]
  "NRXN1", "NRXN2", "NRXN3"    # Neurexins [web:75]
)

# CATEGORY 3: Tissue-restricted ligands — not expressed by T cells
# (Verified against CTCL/Sezary scRNA-seq datasets)
non_tcell_ligands <- c(
  "PRG4",   # Synoviocyte/chondrocyte — absent in T cells
  "C3",     # Liver-secreted complement — T cells do not produce autocrinically
  "TCN1",   # Neutrophil/salivary gland — not in T cells
  "HGF"     # Fibroblast/hepatocyte product — not secreted by T cells [web:114]
)

# CATEGORY 4: Sézary-specific invalid receptors
sezary_receptor_blacklist <- c(
  "DPP4",   # CD26 — absent on Sézary cells (WHO criterion) [web:40]
  "GRM7"    # Glutamate receptor — absent in CD4 T cells [web:66]
)

# CATEGORY 5: Specific LR pair mismatches (wrong receptor class)
# NOTE: CCL4->CCR8 is deliberately NOT here — it is a validated interaction [web:97][web:101]
lr_pair_blacklist <- data.frame(
  ligand   = c("CCL5",  "CXCL10", "CXCL10", "CCL3",  "SPP1",  "SPP1",   "SPP1",   "TNF"),
  receptor = c("GRM7",  "GRM7",   "TLR4",   "CCR4",  "S1PR1", "CCR8",   "PTGER4", "VSIR"),
  stringsAsFactors = FALSE
)

# ─── APPLY ───────────────────────────────────────────────────────────────────
full_ligand_blacklist   <- unique(c(intracellular_blacklist, neuronal_blacklist, non_tcell_ligands))
full_receptor_blacklist <- unique(c(intracellular_blacklist, neuronal_blacklist, sezary_receptor_blacklist))

filtered_data <- filtered_data %>%
  dplyr::filter(
    !ligand.complex   %in% full_ligand_blacklist,
    !receptor.complex %in% full_receptor_blacklist
  ) %>%
  dplyr::anti_join(
    lr_pair_blacklist,
    by = c("ligand.complex" = "ligand", "receptor.complex" = "receptor")
  )

cat("After comprehensive biological filter:", nrow(filtered_data), "\n")

After comprehensive biological filter: 29078 

# ─── AUDIT: Which blacklisted genes were actually present? ────────────────
removed_ligands   <- unique(all_data$ligand.complex[all_data$ligand.complex %in% full_ligand_blacklist])
removed_receptors <- unique(all_data$receptor.complex[all_data$receptor.complex %in% full_receptor_blacklist])

cat("\nBlacklisted ligands found in your data:  ", paste(removed_ligands, collapse=", "), "\n")

Blacklisted ligands found in your data:   PTGS2, TCN1, GNAI2, FADD, NRXN1, NLGN1, PRG4, HGF, C3 

cat("Blacklisted receptors found in your data:", paste(removed_receptors, collapse=", "), "\n")

Blacklisted receptors found in your data: TRAF2, TRADD, DPP4, RIPK1, GRM7, NLGN1, NRXN1 

# ─── APPLY PAIR-SPECIFIC BLACKLIST (run after gene blacklist) ──────────────

lr_pair_blacklist <- data.frame(
  ligand   = c("CXCL10", "CCL3",  "SPP1",  "SPP1",   "SPP1",   "TNF"),
  receptor = c("TLR4",   "CCR4",  "S1PR1", "CCR8",   "PTGER4", "VSIR"),
  stringsAsFactors = FALSE
)

# Check which of these pairs are actually present in your current filtered data
pairs_present <- filtered_data %>%
  dplyr::semi_join(
    lr_pair_blacklist,
    by = c("ligand.complex" = "ligand", "receptor.complex" = "receptor")
  ) %>%
  dplyr::distinct(ligand.complex, receptor.complex)

cat("Pair-specific interactions found and will be removed:\n")

Pair-specific interactions found and will be removed:

print(pairs_present)

[1] ligand.complex   receptor.complex
<0 rows> (or 0-length row.names)

# Apply the anti-join
filtered_data <- filtered_data %>%
  dplyr::anti_join(
    lr_pair_blacklist,
    by = c("ligand.complex" = "ligand", "receptor.complex" = "receptor")
  )

cat("\nFinal interaction count after pair blacklist:", nrow(filtered_data), "\n")

Final interaction count after pair blacklist: 29078 

# ─── APPLY PAIR-SPECIFIC BLACKLIST (run after gene blacklist) ──────────────

# ─── SIGNIFICANCE FILTER ──────────────────────────────────────────────────
significant_interactions <- filtered_data %>%
  dplyr::filter(aggregate_rank <= 0.05)

cat("Significant interactions (rank ≤ 0.05):", nrow(significant_interactions), "\n")

Significant interactions (rank ≤ 0.05): 5644 

# ─── SANITY CHECK: Confirm blacklisted genes are fully gone ───────────────
blacklist_check <- significant_interactions %>%
  dplyr::filter(
    ligand.complex   %in% c(full_ligand_blacklist, "GRM7", "TLR4") |
    receptor.complex %in% c(full_receptor_blacklist, "GRM7", "TLR4")
  )

cat("Blacklisted interactions remaining (must be 0):", nrow(blacklist_check), "\n")

Blacklisted interactions remaining (must be 0): 14 

# ─── DISTRIBUTION CHECK ───────────────────────────────────────────────────
cat("Cell lines represented:", paste(sort(unique(significant_interactions$CellLine)),
                                      collapse = ", "), "\n")

Cell lines represented: L1, L2, L3, L4, L5, L6, L7 

cat("Unique interactions:    ", n_distinct(significant_interactions$Interaction), "\n")

Unique interactions:     327 

cat("Interactions per cell line:\n")

Interactions per cell line:

print(table(significant_interactions$CellLine))

  L1   L2   L3   L4   L5   L6   L7 
 770  153  473 1312 1324  382 1230 

# ─── DIAGNOSE: Find exactly which 14 interactions survived ────────────────

surviving_blacklist <- significant_interactions %>%
  dplyr::filter(
    ligand.complex   %in% c(full_ligand_blacklist, "GRM7", "TLR4") |
    receptor.complex %in% c(full_receptor_blacklist, "GRM7", "TLR4")
  ) %>%
  dplyr::distinct(ligand.complex, receptor.complex, CellLine)

print(surviving_blacklist)

  ligand.complex receptor.complex CellLine
1          ZG16B             TLR4       L5
2         HSPA1A             TLR4       L5
3        HSP90B1             TLR4       L5
4          HMGB1             TLR4       L5
5        HSP90B1             TLR4       L7
6          HMGB1             TLR4       L7
7         HSPA1A             TLR4       L7

# ─── DIAGNOSE: Find exactly which 14 interactions survived ────────────────

# ─── CORRECTED PAIR BLACKLIST: only the specific invalid pair ─────────────
lr_pair_blacklist <- data.frame(
  ligand   = c("CXCL10", "CCL3",  "SPP1",  "SPP1",   "SPP1",   "TNF"),
  receptor = c("TLR4",   "CCR4",  "S1PR1", "CCR8",   "PTGER4", "VSIR"),
  stringsAsFactors = FALSE
)
# This removes CXCL10->TLR4 specifically but KEEPS:
# HMGB1->TLR4   ✅ validated DAMP signaling
# HSPA1A->TLR4  ✅ validated extracellular HSP70
# HSP90B1->TLR4 ✅ validated extracellular HSP90
# ZG16B->TLR4   ✅ keep (OmniPath-supported)

# ─── Also fix the gene-level receptor blacklist: remove TLR4 from it ──────
# TLR4 should NOT be in full_receptor_blacklist — it is a valid surface receptor
full_receptor_blacklist <- full_receptor_blacklist[full_receptor_blacklist != "TLR4"]

# ─── Rebuild the grepl pattern without TLR4 as a receptor ─────────────────
all_blacklisted_genes <- unique(c(
  full_ligand_blacklist,
  full_receptor_blacklist,
  "GRM7"
  # NOTE: TLR4 deliberately excluded — it is a valid DAMP receptor
))

blacklist_pattern <- paste0(
  "\\b(",
  paste(all_blacklisted_genes, collapse = "|"),
  ")\\b"
)

# ─── Reapply clean filter ──────────────────────────────────────────────────
significant_interactions <- filtered_data %>%
  dplyr::filter(aggregate_rank <= 0.05) %>%
  dplyr::filter(
    !grepl(blacklist_pattern, ligand.complex,   perl = TRUE) &
    !grepl(blacklist_pattern, receptor.complex, perl = TRUE)
  ) %>%
  dplyr::anti_join(
    lr_pair_blacklist,
    by = c("ligand.complex" = "ligand", "receptor.complex" = "receptor")
  )

cat("Final significant interactions:", nrow(significant_interactions), "\n")

Final significant interactions: 5644 

# ─── Final verification ────────────────────────────────────────────────────
final_check <- significant_interactions %>%
  dplyr::filter(
    grepl(blacklist_pattern, ligand.complex,   perl = TRUE) |
    grepl(blacklist_pattern, receptor.complex, perl = TRUE)
  )
cat("Remaining blacklisted (must be 0):", nrow(final_check), "\n")

Remaining blacklisted (must be 0): 0 

# Confirm DAMP interactions survived
damp_check <- significant_interactions %>%
  dplyr::filter(ligand.complex %in% c("HMGB1", "HSPA1A", "HSP90B1") &
                receptor.complex == "TLR4") %>%
  dplyr::distinct(ligand.complex, receptor.complex)
cat("DAMP->TLR4 interactions preserved:\n")

DAMP->TLR4 interactions preserved:

print(damp_check)

  ligand.complex receptor.complex
1         HSPA1A             TLR4
2        HSP90B1             TLR4
3          HMGB1             TLR4

# ─── ADD 3 NEW FALSE POSITIVES TO BLACKLIST ───────────────────────────────

intracellular_blacklist <- c(
  intracellular_blacklist,
  "GNAS"    # G-protein alpha-s — cytoplasmic, identical to GNAI2/GNAI3
)

# FABP5->RXRA and NUCB2->ERAP1 are specific pair mismatches
# (FABP5 and NUCB2 can have other valid extracellular roles, so blacklist the pair only)
lr_pair_blacklist <- rbind(
  lr_pair_blacklist,
  data.frame(
    ligand   = c("FABP5", "NUCB2"),
    receptor = c("RXRA",  "ERAP1"),
    stringsAsFactors = FALSE
  )
)

# Rebuild blacklist pattern
all_blacklisted_genes <- unique(c(
  full_ligand_blacklist,
  full_receptor_blacklist,
  "GNAS", "GRM7"
))

blacklist_pattern <- paste0(
  "\\b(",
  paste(all_blacklisted_genes, collapse = "|"),
  ")\\b"
)

# Rerun significance filter with updated lists
significant_interactions <- filtered_data %>%
  dplyr::filter(aggregate_rank <= 0.05) %>%
  dplyr::filter(
    !grepl(blacklist_pattern, ligand.complex,   perl = TRUE) &
    !grepl(blacklist_pattern, receptor.complex, perl = TRUE)
  ) %>%
  dplyr::anti_join(
    lr_pair_blacklist,
    by = c("ligand.complex" = "ligand", "receptor.complex" = "receptor")
  )

cat("Final significant interactions:", nrow(significant_interactions), "\n")

Final significant interactions: 5505 

# Confirm the 3 removed pairs are gone
confirm_removed <- significant_interactions %>%
  dplyr::filter(
    (ligand.complex == "GNAS"  & receptor.complex == "ADCY7") |
    (ligand.complex == "FABP5" & receptor.complex == "RXRA")  |
    (ligand.complex == "NUCB2" & receptor.complex == "ERAP1")
  )
cat("Removed pairs remaining (must be 0):", nrow(confirm_removed), "\n")

Removed pairs remaining (must be 0): 0 

# Confirm VIM, PKM, GPI survived
confirm_kept <- significant_interactions %>%
  dplyr::filter(ligand.complex %in% c("VIM", "PKM", "GPI")) %>%
  dplyr::distinct(ligand.complex, receptor.complex)
cat("Validated moonlighting proteins preserved:\n")

Validated moonlighting proteins preserved:

print(confirm_kept)

  ligand.complex receptor.complex
1            VIM             CD44
2            PKM             CD44
3            GPI             AMFR

---

7 Build Heatmap Score Matrix

# Guard against log10(0) = -Inf with epsilon
heatmap_matrix <- significant_interactions %>%
  dplyr::select(Interaction, CellLine, aggregate_rank) %>%
  dplyr::mutate(score = -log10(aggregate_rank + 1e-10)) %>%
  dplyr::group_by(Interaction, CellLine) %>%
  dplyr::summarise(score = max(score), .groups = "drop") %>%
  tidyr::pivot_wider(
    names_from  = CellLine,
    values_from = score,
    values_fill = 0
  ) %>%
  as.data.frame()

rownames(heatmap_matrix) <- heatmap_matrix$Interaction
heatmap_matrix$Interaction <- NULL
heatmap_matrix <- as.matrix(heatmap_matrix)

cat("Matrix:", nrow(heatmap_matrix), "interactions ×",
    ncol(heatmap_matrix), "cell lines\n")

Matrix: 321 interactions × 7 cell lines

---

8 Filter: Conserved & Strong Interactions

# Conserved: present in ≥5 of 7 cell lines (-log10 score > 1.3 ≡ rank < 0.05)
is_conserved <- rowSums(heatmap_matrix > 1.3) >= 5

# Strong: at least one cell line with rank ≤ 0.001 (score > 3)
is_strong     <- apply(heatmap_matrix, 1, max) > 3

final_matrix  <- heatmap_matrix[is_conserved & is_strong, , drop = FALSE]
final_matrix  <- final_matrix[
  order(rowSums(final_matrix > 1.3), decreasing = TRUE), , drop = FALSE
]

cat("Conserved & strong interactions:", nrow(final_matrix), "\n")

Conserved & strong interactions: 45 

---

9 Conserved Interaction Heatmap

col_fun <- colorRamp2(c(0, 5), c("white", "firebrick"))

ht <- Heatmap(
  as.matrix(final_matrix),
  name             = "-log10(Rank)",
  col              = col_fun,
  cluster_columns  = FALSE,
  cluster_rows     = TRUE,
  show_row_dend    = TRUE,
  show_column_dend = FALSE,
  row_names_gp     = gpar(fontsize = 7),
  column_names_gp  = gpar(fontsize = 10, fontface = "bold"),
  column_title     = "Conserved & Cell Line-Specific Signaling",
  border           = TRUE
)

draw(ht)

pdf("Figure_3_15_Conserved_Heatmap.pdf", width = 7, height = 7)

  draw(ht)
dev.off()

png 
  2 

png("Figure_3_15_Conserved_Heatmap.png", width = 1800, height = 2200, res = 300)
  draw(ht)

dev.off()

png 
  2 

cat("✓ Conserved heatmap saved\n")

✓ Conserved heatmap saved

---

10 Top 50 Strongest Interactions Heatmap

max_scores    <- apply(heatmap_matrix, 1, max)
n_top         <- min(50, nrow(heatmap_matrix))
top_50_names  <- names(sort(max_scores, decreasing = TRUE))[seq_len(n_top)]
top50_matrix  <- as.matrix(heatmap_matrix[top_50_names, , drop = FALSE])
top50_matrix  <- top50_matrix[
  order(rowSums(top50_matrix > 1.3), decreasing = TRUE), , drop = FALSE
]

ht_top50 <- Heatmap(
  top50_matrix,
  name             = "-log10(Rank)",
  col              = col_fun,
  cluster_columns  = FALSE,
  cluster_rows     = TRUE,
  show_row_dend    = TRUE,
  show_column_dend = FALSE,
  row_names_gp     = gpar(fontsize = 8),
  column_names_gp  = gpar(fontsize = 10, fontface = "bold"),
  column_title     = "Top 50 Strongest Interactions",
  border           = TRUE
)

draw(ht_top50)

pdf("Figure_Top50_Strongest_Interactions.pdf", width = 8, height = 10)

  draw(ht_top50)
dev.off()

png 
  2 

png("Figure_Top50_Strongest_Interactions.png", width = 2400, height = 3000, res = 300)
  draw(ht_top50)

dev.off()

png 
  2 

cat("✓ Top 50 heatmap saved\n")

✓ Top 50 heatmap saved

---

11 Master Ligands

ligand_activity <- filtered_data %>%
  dplyr::group_by(ligand.complex) %>%
  dplyr::summarise(
    interactions = n(),
    receptors    = n_distinct(receptor.complex),
    cell_lines   = n_distinct(CellLine),
    .groups      = "drop"
  ) %>%
  dplyr::arrange(desc(interactions))

print(head(ligand_activity, 20))

# A tibble: 20 × 4
   ligand.complex interactions receptors cell_lines
   <chr>                 <int>     <int>      <int>
 1 CALM1                  1348         8          7
 2 B2M                    1184         9          7
 3 HLA-A                  1101         8          7
 4 TGFB1                  1003         7          7
 5 HLA-B                   951         6          7
 6 LRPAP1                  929         6          7
 7 FN1                     798        10          5
 8 HLA-C                   712         6          7
 9 PSEN1                   707         4          7
10 ICAM1                   696         4          7
11 LTA                     638         3          7
12 APP                     488         2          7
13 TNFSF12                 473         3          7
14 TNF                     463         3          6
15 SELPLG                  447         3          7
16 PTDSS1                  435         2          7
17 ADAM17                  414         4          7
18 TNFSF10                 399         3          7
19 CCL5                    395         4          7
20 GNAS                    358         4          7

---

12 Master Ligands Plot

ggplot(head(ligand_activity, 20),
       aes(x = reorder(ligand.complex, interactions), y = interactions)) +
  geom_bar(stat = "identity", fill = "firebrick") +
  coord_flip() +
  theme_classic(base_size = 12) +
  labs(x = "Ligand", y = "Number of Interactions",
       title = "Master Signaling Ligands — Sézary Syndrome Cell Lines")

ggsave("Figure_Master_Ligands.pdf", width = 8, height = 6)
ggsave("Figure_Master_Ligands.png", width = 8, height = 6, dpi = 300)

---

13 Conserved Interactions Table

conserved_threshold <- 4

conserved_interactions <- as.data.frame(final_matrix) %>%
  dplyr::mutate(
    Interaction  = rownames(.),
    NumCellLines = rowSums(dplyr::across(where(is.numeric)) > 1.3)
  ) %>%
  dplyr::filter(NumCellLines >= conserved_threshold) %>%
  dplyr::arrange(desc(NumCellLines)) %>%
  dplyr::select(Interaction, NumCellLines, dplyr::everything())

cat("Conserved interactions (≥", conserved_threshold, "cell lines):",
    nrow(conserved_interactions), "\n")

Conserved interactions (≥ 4 cell lines): 45 

print(head(conserved_interactions[, c("Interaction", "NumCellLines")], 20))

                  Interaction NumCellLines
APP -> CD74       APP -> CD74            7
B2M -> CD3D       B2M -> CD3D            7
B2M -> CD3G       B2M -> CD3G            7
B2M -> TFRC       B2M -> TFRC            7
CALM1 -> HMMR   CALM1 -> HMMR            7
CALM1 -> KCNQ5 CALM1 -> KCNQ5            7
CALM1 -> PTPRA CALM1 -> PTPRA            7
CALM3 -> KCNQ5 CALM3 -> KCNQ5            7
CD48 -> CD2       CD48 -> CD2            7
CD58 -> CD2       CD58 -> CD2            7
HLA-A -> APLP2 HLA-A -> APLP2            7
HLA-A -> CD3D   HLA-A -> CD3D            7
HLA-A -> CD3G   HLA-A -> CD3G            7
HLA-B -> CANX   HLA-B -> CANX            7
HLA-B -> CD3D   HLA-B -> CD3D            7
HLA-C -> CD3D   HLA-C -> CD3D            7
HSPA8 -> LDLR   HSPA8 -> LDLR            7
PKM -> CD44       PKM -> CD44            7
B2M -> CD247     B2M -> CD247            6
CALM1 -> KCNN4 CALM1 -> KCNN4            6

write.csv(conserved_interactions, "Conserved_Interactions_SS.csv", row.names = FALSE)

---

14 Cell Line–Specific Interactions

specific_interactions <- as.data.frame(final_matrix) %>%
  dplyr::mutate(
    Interaction  = rownames(.),
    NumCellLines = rowSums(dplyr::across(where(is.numeric)) > 1.3)
  ) %>%
  dplyr::filter(NumCellLines <= 5) %>%
  dplyr::arrange(NumCellLines) %>%
  dplyr::select(Interaction, NumCellLines, dplyr::everything())

cat("Cell line-specific interactions:", nrow(specific_interactions), "\n")

Cell line-specific interactions: 11 

print(head(specific_interactions[, c("Interaction", "NumCellLines")], 20))

                              Interaction NumCellLines
CCL5 -> SDC4                 CCL5 -> SDC4            5
CXCL10 -> SDC4             CXCL10 -> SDC4            5
FN1 -> CD44                   FN1 -> CD44            5
FN1 -> ITGA4_ITGB1     FN1 -> ITGA4_ITGB1            5
FN1 -> ITGA4_ITGB7     FN1 -> ITGA4_ITGB7            5
FN1 -> ITGAV_ITGB1     FN1 -> ITGAV_ITGB1            5
HLA-A -> KIR3DL2         HLA-A -> KIR3DL2            5
HLA-DRA -> CD4             HLA-DRA -> CD4            5
LTA -> TNFRSF14           LTA -> TNFRSF14            5
LTB -> TNFRSF1A           LTB -> TNFRSF1A            5
TNFSF10 -> TNFRSF10B TNFSF10 -> TNFRSF10B            5

write.csv(specific_interactions, "CellLineSpecific_Interactions.csv", row.names = FALSE)

---

15 UpSetR — Interaction Sharing

upset_df        <- as.data.frame(heatmap_matrix > 1.3) * 1L
available_lines <- intersect(
  c("L1","L2","L3","L4","L5","L6","L7"),
  colnames(upset_df)
)

upset(
  upset_df,
  sets              = rev(available_lines),
  order.by          = "freq",
  keep.order        = TRUE,
  sets.bar.color    = "firebrick",
  main.bar.color    = "black",
  matrix.color      = "black",
  mainbar.y.label   = "Number of Shared Interactions",
  sets.x.label      = "Total Significant Interactions",
  text.scale        = c(1.5, 1.2, 1.2, 1.2, 1.5, 1.2),
  mb.ratio          = c(0.6, 0.4)
)

pdf("Figure_3_15B_UpSetR.pdf", width = 10, height = 7, onefile = FALSE)
  upset(upset_df, sets = rev(available_lines), order.by = "freq",
        keep.order = TRUE, sets.bar.color = "firebrick",
        mainbar.y.label = "Number of Shared Interactions",
        sets.x.label    = "Total Significant Interactions",
        text.scale = c(1.5, 1.2, 1.2, 1.2, 1.5, 1.2))

dev.off()

png 
  2 

png("Figure_3_15B_UpSetR.png", width = 3000, height = 2100, res = 300)

  upset(upset_df, sets = rev(available_lines), order.by = "freq",
        keep.order = TRUE, sets.bar.color = "firebrick",
        mainbar.y.label = "Number of Shared Interactions",
        sets.x.label    = "Total Significant Interactions",
        text.scale = c(1.5, 1.2, 1.2, 1.2, 1.5, 1.2))
dev.off()

png 
  2 

cat("✓ UpSetR plot saved\n")

✓ UpSetR plot saved

---

16 Export & Summary

write.csv(ligand_activity, "Master_Ligands_SS.csv", row.names = FALSE)

cat("=== ANALYSIS SUMMARY ===\n")

=== ANALYSIS SUMMARY ===

cat("Total LIANA records:              ", nrow(all_data),               "\n")

Total LIANA records:               48958 

cat("After OmniPath filter:            ", nrow(filtered_data),          "\n")

After OmniPath filter:             29078 

cat("Significant (rank ≤ 0.05):        ", nrow(significant_interactions),"\n")

Significant (rank ≤ 0.05):         5505 

cat("Conserved & strong (final):       ", nrow(final_matrix),           "\n")

Conserved & strong (final):        45 

cat("Conserved (≥4 lines):             ", nrow(conserved_interactions),  "\n")

Conserved (≥4 lines):              45 

cat("Cell line-specific (≤2 lines):    ", nrow(specific_interactions),   "\n")

Cell line-specific (≤2 lines):     11 

---

17 Session Info

sessionInfo()

R version 4.5.2 (2025-10-31)
Platform: x86_64-pc-linux-gnu
Running under: Ubuntu 24.04.3 LTS

Matrix products: default
BLAS:   /usr/lib/x86_64-linux-gnu/blas/libblas.so.3.12.0 
LAPACK: /usr/lib/x86_64-linux-gnu/lapack/liblapack.so.3.12.0  LAPACK version 3.12.0

locale:
 [1] LC_CTYPE=en_GB.UTF-8       LC_NUMERIC=C               LC_TIME=fr_FR.UTF-8        LC_COLLATE=en_GB.UTF-8    
 [5] LC_MONETARY=fr_FR.UTF-8    LC_MESSAGES=en_GB.UTF-8    LC_PAPER=fr_FR.UTF-8       LC_NAME=C                 
 [9] LC_ADDRESS=C               LC_TELEPHONE=C             LC_MEASUREMENT=fr_FR.UTF-8 LC_IDENTIFICATION=C       

time zone: Europe/Paris
tzcode source: system (glibc)

attached base packages:
[1] grid      stats     graphics  grDevices utils     datasets  methods   base     

other attached packages:
[1] liana_0.1.14          stringr_1.6.0         UpSetR_1.4.0          circlize_0.4.17       ComplexHeatmap_2.26.1
[6] ggplot2_4.0.2         tidyr_1.3.2           dplyr_1.2.0          

loaded via a namespace (and not attached):
  [1] RColorBrewer_1.1-3          rstudioapi_0.18.0           jsonlite_2.0.0              shape_1.4.6.1              
  [5] magrittr_2.0.4              magick_2.9.0                farver_2.1.2                rmarkdown_2.30             
  [9] ragg_1.5.0                  GlobalOptions_0.1.3         fs_1.6.6                    vctrs_0.7.1                
 [13] memoise_2.0.1               Cairo_1.7-0                 htmltools_0.5.9             S4Arrays_1.10.1            
 [17] progress_1.2.3              curl_7.0.0                  BiocNeighbors_2.4.0         cellranger_1.1.0           
 [21] SparseArray_1.10.8          parallelly_1.46.1           basilisk_1.23.0             plyr_1.8.9                 
 [25] httr2_1.2.2                 lubridate_1.9.5             cachem_1.1.0                igraph_2.2.2               
 [29] lifecycle_1.0.5             iterators_1.0.14            pkgconfig_2.0.3             rsvd_1.0.5                 
 [33] Matrix_1.7-4                R6_2.6.1                    fastmap_1.2.0               MatrixGenerics_1.22.0      
 [37] future_1.69.0               clue_0.3-67                 digest_0.6.39               colorspace_2.1-2           
 [41] S4Vectors_0.48.0            dqrng_0.4.1                 irlba_2.3.7                 textshaping_1.0.4          
 [45] GenomicRanges_1.62.1        RSQLite_2.4.6               beachmat_2.26.0             labeling_0.4.3             
 [49] filelock_1.0.3              progressr_0.18.0            timechange_0.4.0            httr_1.4.8                 
 [53] abind_1.4-8                 compiler_4.5.2              bit64_4.6.0-1               withr_3.0.2                
 [57] doParallel_1.0.17           S7_0.2.1                    backports_1.5.0             BiocParallel_1.44.0        
 [61] DBI_1.3.0                   logger_0.4.1                OmnipathR_3.19.1            R.utils_2.13.0             
 [65] rappdirs_0.3.4              DelayedArray_0.36.0         sessioninfo_1.2.3           bluster_1.20.0             
 [69] rjson_0.2.23                tools_4.5.2                 otel_0.2.0                  zip_2.3.3                  
 [73] future.apply_1.20.2         R.oo_1.27.1                 glue_1.8.0                  checkmate_2.3.4            
 [77] cluster_2.1.8.2             generics_0.1.4              gtable_0.3.6                tzdb_0.5.0                 
 [81] R.methodsS3_1.8.2           hms_1.1.4                   utf8_1.2.6                  metapod_1.18.0             
 [85] ScaledMatrix_1.18.0         BiocSingular_1.26.1         sp_2.2-1                    xml2_1.5.2                 
 [89] XVector_0.50.0              BiocGenerics_0.56.0         foreach_1.5.2               pillar_1.11.1              
 [93] limma_3.66.0                spam_2.11-3                 later_1.4.7                 lattice_0.22-9             
 [97] bit_4.6.0                   tidyselect_1.2.1            locfit_1.5-9.12             SingleCellExperiment_1.32.0
[101] scuttle_1.20.0              knitr_1.51                  gridExtra_2.3               IRanges_2.44.0             
[105] Seqinfo_1.0.0               edgeR_4.8.2                 SummarizedExperiment_1.40.0 stats4_4.5.2               
[109] xfun_0.56                   Biobase_2.70.0              statmod_1.5.1               matrixStats_1.5.0          
[113] stringi_1.8.7               yaml_2.3.12                 evaluate_1.0.5              codetools_0.2-20           
[117] tcltk_4.5.2                 tibble_3.3.1                cli_3.6.5                   systemfonts_1.3.1          
[121] reticulate_1.45.0           dichromat_2.0-0.1           Rcpp_1.1.1                  readxl_1.4.5               
[125] globals_0.19.0              dir.expiry_1.18.0           png_0.1-8                   XML_3.99-0.22              
[129] parallel_4.5.2              readr_2.2.0                 blob_1.3.0                  prettyunits_1.2.0          
[133] dotCall64_1.2               scran_1.38.0                listenv_0.10.0              scales_1.4.0               
[137] SeuratObject_5.3.0          purrr_1.2.1                 crayon_1.5.3                GetoptLong_1.1.0           
[141] rlang_1.1.7                 rvest_1.0.5