Liana analysis L1-Thesis

1 load libraries

2 load seurat object

# 1. Reload clean object
L1 <- readRDS("../../../0-RDS_Cell_lines/L1_clustered.rds")
DefaultAssay(L1) <- "SCT"

Idents(L1) <- "SCT_snn_res.0.2"

3 Error Fix

my_GetAssayData <- function(
  object,
  assay = NULL,
  layer = NULL,
  slot = NULL,
  ...
) {
  # Map old 'slot' to new 'layer'
  if (is.null(layer) && !is.null(slot)) {
    layer <- slot
  }
  if (is.null(layer)) {
    layer <- "data"
  }

  # If a Seurat object is provided, get the assay object
  if (inherits(object, "Seurat")) {
    if (is.null(assay)) assay <- Seurat::DefaultAssay(object)
    assay_obj <- object[[assay]]
  } else {
    assay_obj <- object
  }

  # Seurat v5 assays: use LayerData
  if (inherits(assay_obj, "Assay5") || inherits(assay_obj, "StdAssay")) {
    return(SeuratObject::LayerData(assay_obj, layer = layer))
  }

  # Seurat v4 assays: fall back to slots
  if (inherits(assay_obj, "Assay")) {
    return(methods::slot(assay_obj, layer))
  }

  stop("Unsupported object class: ", paste(class(object), collapse = ", "))
}

assignInNamespace("GetAssayData", my_GetAssayData, ns = "SeuratObject")

4 Liana analysis using scPubr

# Run LIANA with multiple methods
liana_output <- liana::liana_wrap(sce =  L1,
                                  method = c("natmi", "connectome", "logfc", "sca", "cellphonedb", "CellChat"),
                                  idents_col = NULL,
                                  verbose = FALSE,
                                  assay = "SCT")

# Compute consensus aggregate ranking
liana_aggregate <- liana::liana_aggregate(liana_output)

4.1 Save Results

# Create output directory
dir.create("LIANA_L1", showWarnings = FALSE)

# Save R object (for future use)
saveRDS(liana_output, "LIANA_L1/liana_L1_output.rds")

# Save aggregate results (for comparative analysis across all cell lines)
write.csv(liana_aggregate, 
          "LIANA_L1/liana_L1_aggregate_results.csv", 
          row.names = FALSE)

# Save individual methods
for (name in names(liana_output)) {
    write.csv(liana_output[[name]], 
              paste0("LIANA_L1/liana_L1_", name, "_results.csv"), 
              row.names = FALSE)
}

cat("✓ Results saved to LIANA_L1/\n")

✓ Results saved to LIANA_L1/

4.2 Supplementary Figure: Detailed Dotplot

library(ggplot2)

# Top 10 interactions by aggregate_rank
p_dotplot <- SCpubr::do_LigandReceptorPlot(
    liana_output = liana_output,
    arrange_interactions_by = "aggregate_rank",
    top_interactions = 10
)

p_dotplot

# Save as PDF
ggsave(
  filename = "LIANA_L1/L1_Dotplot_Top10_AggregateRank.pdf",
  plot = p_dotplot,
  width = 10,
  height = 8
)

# Save as PNG (set dpi explicitly — otherwise figures look amateurish)
ggsave(
  filename = "LIANA_L1/L1_Dotplot_Top10_AggregateRank.png",
  plot = p_dotplot,
  width = 10,
  height = 8,
  dpi = 300
)

cat("✓ Dotplot saved as PDF and PNG\n")

✓ Dotplot saved as PDF and PNG

5 Main Figure: Chord Diagrams

# Generate chord diagrams (TOP 25 interactions by aggregate_rank)
chord_output <- SCpubr::do_LigandReceptorPlot(
    liana_output = liana_output,
    top_interactions = 10,
    arrange_interactions_by = "aggregate_rank",
    compute_ChordDiagrams = TRUE
)

# Display chord diagrams
chord_output$chord_total_interactions
chord_output$chord_ligand_receptor

5.1 Save Chord Diagrams

# Save Chord Diagram 1: Total interactions by cluster
pdf("LIANA_L1/L1_Chord_Total_Interactions.pdf", width = 10, height = 8)
print(chord_output$chord_total_interactions)
dev.off()

null device 
          1 

# Save Chord Diagram 2: Specific ligand-receptor pairs
pdf("LIANA_L1/L1_Chord_Ligand_Receptor.pdf", width = 10, height = 8)
print(chord_output$chord_ligand_receptor)
dev.off()

null device 
          1 

# Save Chord Diagram 1: Total interactions by cluster
pdf("LIANA_L1/L1_Chord_Total_Interactions.png", width = 10, height = 8)
print(chord_output$chord_total_interactions)
dev.off()

null device 
          1 

# Save Chord Diagram 2: Specific ligand-receptor pairs
pdf("LIANA_L1/L1_Chord_Ligand_Receptor.png", width = 10, height = 8)
print(chord_output$chord_ligand_receptor)
dev.off()

null device 
          1 

cat("✓ Chord diagrams saved\n")

✓ Chord diagrams saved

6 Summary Table

# Top 10 interactions for manuscript table
top10_table <- liana_aggregate %>%
  arrange(aggregate_rank) %>%
  select(source, target, ligand.complex, receptor.complex, 
         aggregate_rank, natmi.edge_specificity, sca.LRscore) %>%
  head(10)

# Display
print(top10_table)

# A tibble: 10 × 7
   source target ligand.complex receptor.complex aggregate_rank natmi.edge_specificity sca.LRscore
   <chr>  <chr>  <chr>          <chr>                     <dbl>                  <dbl>       <dbl>
 1 5      5      TNF            ICOS                0.000000518                 0.122        0.792
 2 5      1      PTGS2          PTGDR2              0.00000121                  0.241        0.675
 3 5      3      PTGS2          PTGDR2              0.0000109                   0.212        0.660
 4 5      4      TNF            VSIR                0.0000234                   0.0832       0.766
 5 4      4      TCN1           CD320               0.0000263                   0.233        0.699
 6 5      3      TNF            TRAF2               0.0000263                   0.0824       0.698
 7 5      0      PTGS2          PTGDR2              0.0000328                   0.144        0.616
 8 4      2      TCN1           CD320               0.0000364                   0.190        0.677
 9 4      4      CD99           CD81                0.0000403                   0.0931       0.773
10 5      3      TNF            SEMA4C              0.0000403                   0.0804       0.664

# Save as CSV
write.csv(top10_table, 
          "LIANA_L1/L1_Top10_Interactions_Table.csv", 
          row.names = FALSE)

cat("✓ Summary table saved\n")

✓ Summary table saved

# Print statistics
cat("\n=== SUMMARY STATISTICS ===\n")

=== SUMMARY STATISTICS ===

cat("Total interactions with aggregate_rank ≤ 0.05:", 
    sum(liana_aggregate$aggregate_rank <= 0.05), "\n")

Total interactions with aggregate_rank ≤ 0.05: 1241 

cat("Total interactions with aggregate_rank ≤ 0.01:", 
    sum(liana_aggregate$aggregate_rank <= 0.01), "\n")

Total interactions with aggregate_rank ≤ 0.01: 442 

6.1 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] stats     graphics  grDevices utils     datasets  methods   base     

other attached packages:
[1] ggplot2_4.0.2      SCpubr_1.1.1.9000  liana_0.1.14       dplyr_1.2.0        Seurat_5.4.0       SeuratObject_5.3.0
[7] sp_2.2-0          

loaded via a namespace (and not attached):
  [1] fs_1.6.6                    matrixStats_1.5.0           spatstat.sparse_3.1-0       lubridate_1.9.5            
  [5] httr_1.4.7                  RColorBrewer_1.1-3          doParallel_1.0.17           tools_4.5.2                
  [9] sctransform_0.4.3           backports_1.5.0             utf8_1.2.6                  R6_2.6.1                   
 [13] lazyeval_0.2.2              uwot_0.2.4                  GetoptLong_1.1.0            withr_3.0.2                
 [17] prettyunits_1.2.0           gridExtra_2.3               progressr_0.18.0            textshaping_1.0.4          
 [21] cli_3.6.5                   Biobase_2.70.0              spatstat.explore_3.7-0      fastDummies_1.7.5          
 [25] labeling_0.4.3              sass_0.4.10                 S7_0.2.1                    spatstat.data_3.1-9        
 [29] readr_2.1.6                 ggridges_0.5.7              pbapply_1.7-4               systemfonts_1.3.1          
 [33] R.utils_2.13.0              scater_1.38.0               dichromat_2.0-0.1           parallelly_1.46.1          
 [37] sessioninfo_1.2.3           limma_3.66.0                readxl_1.4.5                rstudioapi_0.18.0          
 [41] RSQLite_2.4.6               generics_0.1.4              shape_1.4.6.1               ica_1.0-3                  
 [45] spatstat.random_3.4-4       zip_2.3.3                   Matrix_1.7-4                ggbeeswarm_0.7.3           
 [49] S4Vectors_0.48.0            logger_0.4.1                abind_1.4-8                 R.methodsS3_1.8.2          
 [53] lifecycle_1.0.5             yaml_2.3.12                 edgeR_4.8.2                 SummarizedExperiment_1.40.0
 [57] SparseArray_1.10.8          Rtsne_0.17                  grid_4.5.2                  blob_1.3.0                 
 [61] promises_1.5.0              dqrng_0.4.1                 crayon_1.5.3                dir.expiry_1.18.0          
 [65] miniUI_0.1.2                lattice_0.22-7              beachmat_2.26.0             cowplot_1.2.0              
 [69] pillar_1.11.1               knitr_1.51                  ComplexHeatmap_2.26.1       metapod_1.18.0             
 [73] GenomicRanges_1.62.1        tcltk_4.5.2                 rjson_0.2.23                future.apply_1.20.1        
 [77] codetools_0.2-20            glue_1.8.0                  spatstat.univar_3.1-6       data.table_1.18.2.1        
 [81] vctrs_0.7.1                 png_0.1-8                   spam_2.11-3                 cellranger_1.1.0           
 [85] gtable_0.3.6                assertthat_0.2.1            cachem_1.1.0                OmnipathR_3.19.1           
 [89] xfun_0.56                   S4Arrays_1.10.1             mime_0.13                   Seqinfo_1.0.0              
 [93] survival_3.8-3              SingleCellExperiment_1.32.0 iterators_1.0.14            statmod_1.5.1              
 [97] bluster_1.20.0              fitdistrplus_1.2-6          ROCR_1.0-12                 nlme_3.1-168               
[101] bit64_4.6.0-1               progress_1.2.3              filelock_1.0.3              RcppAnnoy_0.0.23           
[105] bslib_0.10.0                irlba_2.3.7                 vipor_0.4.7                 KernSmooth_2.23-26         
[109] otel_0.2.0                  colorspace_2.1-2            BiocGenerics_0.56.0         DBI_1.2.3                  
[113] tidyselect_1.2.1            bit_4.6.0                   compiler_4.5.2              curl_7.0.0                 
[117] rvest_1.0.5                 httr2_1.2.2                 BiocNeighbors_2.4.0         xml2_1.5.2                 
[121] DelayedArray_0.36.0         plotly_4.12.0               checkmate_2.3.4             scales_1.4.0               
[125] lmtest_0.9-40               rappdirs_0.3.4              stringr_1.6.0               digest_0.6.39              
[129] goftest_1.2-3               spatstat.utils_3.2-1        rmarkdown_2.30              basilisk_1.23.0            
[133] XVector_0.50.0              htmltools_0.5.9             pkgconfig_2.0.3             sparseMatrixStats_1.22.0   
[137] MatrixGenerics_1.22.0       fastmap_1.2.0               rlang_1.1.7                 GlobalOptions_0.1.3        
[141] htmlwidgets_1.6.4           shiny_1.12.1                jquerylib_0.1.4             farver_2.1.2               
[145] zoo_1.8-15                  jsonlite_2.0.0              BiocParallel_1.44.0         R.oo_1.27.1                
[149] BiocSingular_1.26.1         magrittr_2.0.4              scuttle_1.20.0              dotCall64_1.2              
[153] patchwork_1.3.2             Rcpp_1.1.1                  viridis_0.6.5               reticulate_1.44.1          
[157] stringi_1.8.7               MASS_7.3-65                 plyr_1.8.9                  parallel_4.5.2             
[161] listenv_0.10.0              ggrepel_0.9.6               forcats_1.0.1               deldir_2.0-4               
[165] splines_4.5.2               tensor_1.5.1                hms_1.1.4                   circlize_0.4.17            
[169] locfit_1.5-9.12             igraph_2.2.1                spatstat.geom_3.7-0         RcppHNSW_0.6.0             
[173] reshape2_1.4.5              stats4_4.5.2                ScaledMatrix_1.18.0         XML_3.99-0.20              
[177] evaluate_1.0.5              scran_1.38.0                tzdb_0.5.0                  foreach_1.5.2              
[181] httpuv_1.6.16               RANN_2.6.2                  tidyr_1.3.2                 purrr_1.2.1                
[185] polyclip_1.10-7             future_1.69.0               clue_0.3-66                 scattermore_1.2            
[189] rsvd_1.0.5                  xtable_1.8-4                RSpectra_0.16-2             later_1.4.5                
[193] ragg_1.5.0                  viridisLite_0.4.3           tibble_3.3.1                beeswarm_0.4.0             
[197] memoise_2.0.1               IRanges_2.44.0              cluster_2.1.8.1             timechange_0.4.0           
[201] globals_0.19.0