Run liana for Ligand Receptor Interaction Prediction

Run LIANA for Ligand-Receptor Interaction Prediction in R

As input for CrossTalkeR we need tables with predicted ligand-receptor interactions for different sample conditions. We recommend using the liana framework, which provides multiple ligand-receptor interaction inference methods and databases. The analysis can be performed on Seurat or SingleCellExperiment objects. Here, we show how to run liana with an Seurat object as input, using CellphoneDB as method and the Consensus database:

library(Seurat)
library(dplyr)
library(liana)

seurat_object <- readRDS("/path/to/Seurat/object.Rds")

outpath <- "/path/to/save/results/"
pval_threshold <- 0.05

seurat_object_list <- SplitObject(seurat_object, split.by = "condition")
for (condition in names(seurat_object_list)) {
  sub_object <- seurat_object_list[[condition]]

  liana_results <- liana_wrap(sub_object,
    method = "cellphonedb",
    resource = c("Consensus"),
    expr_prop = 0.1
  )

  liana_results <- liana_results %>%
    filter(pvalue < pval_threshold) %>%
    select(source, target, ligand, receptor.complex, lr.mean) %>%
    rename(gene_A = ligand, gene_B = receptor.complex, MeanLR = lr.mean)
  liana_results$type_gene_A = 'Ligand'
  liana_results$type_gene_B = 'Receptor'

  write.csv(liana_results, paste0(outpath, condition, "_LR.csv"))
}

First, we are importing all necessary libraries and loading the Seurat object. Since we want to predict ligand-receptor interactions for each sample condition in our object, we are splitting the object by using the corresponding meta data column. Then we run liana for each of the condition, and filter the results by the p-value provided. Before saving the result table, we reformat the results to match CrossTalkeRs input format by renaming and adding gene type information columns. The resulting tables can be directly used for a CrossTalkeR analysis.

Run LIANA for Ligand-Receptor Interaction Prediction in Python

We can also use the python version of the package liana. This package provides different methods and databases for the analysis of intercellular interactions based on scRNA-seq data in the AnnData format. In this Tutorial, we describe our recommendation for using the liana package in python.

Import Necessary libraries

First of all, the libraries necessary for using liana must be loaded:

import scanpy as sc
import liana as li
import pandas as pd
import os

We assume here that the necessary libraries are already installed.

Load scRNA-seq Data

The last step before perfoming the analysis is to load the AnnData object with the scRNA-seq data

data = sc.read_h5ad("/path/to/AnnData/file.h5ad")
data.raw = data

Execution of Ligand-Receptor Interaction Analysis

As mentioned at the beginning, there are different methods and databases implemented in the liana package. We recommend using the CellphoneDB analysis in combination with the consensus database with intercellular interactions. Note that metadata fields for the available sample/conditions (field_with_condition) and for the cell annotations (field_with_cell_annotation) are still specified by the user.

for i in set(data.obs['field_with_condition']):
    print(i)
    lr=li.method.cellphonedb(data[data.obs['field_with_condition']==i],
                          groupby='field_with_cell_annotation',
                          expr_prop=0.1,
                          verbose=True,
                          resource_name='consensus',
                          inplace=False)
    lr.to_csv(f"{i}_lr_liana_consensus_unfiltered.csv")

Filtering and Formatting of the Results

Finally, the result tables are filtered by the CellphoneDB p-value (here in the example p-val < 0.05) and converted into the CrossTalkeR input format. By default, the interaction table used for CrossTalkeR should contain the columns “source”, “target”, “type_gene_A”, “gene_A”, “type_gene_B”, “gene_B” and “MeanLR”.

data = {}
for i in os.listdir():
    if i.endswith('lr_liana_consensus_unfiltered.csv'):
        evfull = pd.read_csv(i)
        evfull = evfull.loc[:,['ligand','receptor','source','target','lr_means','cellphone_pvals']]
        evfull['type_gene_A'] = 'Ligand'
        evfull['type_gene_B'] = 'Receptor'
        evfull['gene_A'] = evfull['ligand']
        evfull['gene_B'] = evfull['receptor']
        evfull['MeanLR'] = evfull['lr_means']
        evfull.loc[list(evfull.cellphone_pvals.to_numpy()<=0.05),:].to_csv(f'{k}_lr_ready.csv')
        evfull = evfull.loc[:, ['source', 'target', 'type_gene_A', 'type_gene_B', 'gene_A', 'gene_B', 'MeanLR']]
        k=i[0:i.find('_lr_')]
        data[k]=os.path.abspath(f'{k}_lr_ready.csv')