saeyslab / nichenetr

NicheNet: predict active ligand-target links between interacting cells
469 stars 117 forks source link

Function generate_info_tables return an error #287

Closed pariaaliour closed 2 months ago

pariaaliour commented 2 months ago

Dear Saeyslab team, I am running nichener for my single cell analysis. I would like to prioritize the ligand-receptor pairs. However, I am getting an error. as below. Not sure if it's a bug or I am doing sth by mistake.

seuratObj <- readRDS(seuratObj.rds")
Idents(seuratObj) <- "cluster_id"
DefaultAssay(seuratObj) <- "RNA"
seuratObj <- alias_to_symbol_seurat(seuratObj, "human")
table(seuratObj$cluster_id)

print("Read in NicheNet’s networks")
organism = "human"
lr_network = readRDS("/home/paria/scratch/als-project/analysis/ocu_med_sc/rpca/nichenetr/lr_network_human_21122021.rds")
ligand_target_matrix = readRDS("/home/paria/scratch/als-project/analysis/ocu_med_sc/rpca/nichenetr/ligand_target_matrix_nsga2r_final.rds")
weighted_networks = readRDS("/home/paria/scratch/als-project/analysis/ocu_med_sc/rpca/nichenetr/weighted_networks_nsga2r_final.rds")
weighted_networks_lr = weighted_networks$lr_sig %>% inner_join(lr_network, by = c("from","to"))

# receiver
receiver="Mic1"
expressed_genes_receiver = get_expressed_genes(receiver, seuratObj, pct = 0.05)
all_receptors <- unique(lr_network$to)  
expressed_receptors = intersect(all_receptors,expressed_genes_receiver)

## sender
sender_celltypes = c("ExInVentral1", "Mic1", "Olig3", "Olig5", "GMAstrocyte", "WMAstrocyte", "OPC1", "OPC2", "EC", "MN")
potential_ligands <- lr_network %>% filter(to %in% expressed_receptors) %>% pull(from) %>% unique()
list_expressed_genes_sender = sender_celltypes %>% unique() %>% lapply(get_expressed_genes, seuratObj, 0.05)
# lapply to get the expressed genes of every sender cell type separately here
expressed_genes_sender = list_expressed_genes_sender %>% unlist() %>% unique()
length(expressed_genes_sender)

#Define a gene set of interest
print("Define a gene set of interest")
seurat_obj_receiver= subset(seuratObj, idents = receiver)

condition_oi = "als"
condition_reference = "con"
DE_table_receiver = FindMarkers(object = seurat_obj_receiver, ident.1 = condition_oi, ident.2 = condition_reference, group.by = "group_id", min.pct = 0.05) %>% rownames_to_column("gene")
geneset_oi = DE_table_receiver %>% dplyr::filter(p_val_adj <= 0.05 & abs(avg_log2FC) >= 0.25) %>% pull(gene)
geneset_oi = geneset_oi %>% .[. %in% rownames(ligand_target_matrix)]
length(geneset_oi)

all_ligands <- unique(lr_network$from)
all_receptors <- unique(lr_network$to)

expressed_ligands <- intersect(all_ligands, expressed_genes_sender)
expressed_receptors <- intersect(all_receptors, expressed_genes_receiver)

lr_network_filtered <-  lr_network %>% filter(from %in% expressed_ligands & to %in% expressed_receptors)

info_tables <- generate_info_tables(seuratObj,
                                      celltype_colname = "cluster_id",
                                      senders_oi = sender_celltypes,
                                      receivers_oi = receiver,
                                      lr_network = lr_network_filtered,
                                      condition_colname = "group_id",
                                      condition_oi = condition_oi,
                                      condition_reference = condition_reference,
                                      scenario = "one_condition")
Error in `rename()`:
! Can't rename columns that don't exist.
✖ Column `cluster` doesn't exist.
Run `rlang::last_trace()` to see where the error occurred.
Warning message:
In generate_info_tables(seuratObj, celltype_colname = "cluster_id",  :
  condition_* arguments are provided but the 'one_condition' scenario is selected. Only cells from the condition_oi will be used to calculate cell type specificity, condition specificity will not be calculated.

Cluster_id and group_id cluster are correctly defined. Not sure why it returns an error.

 table(seuratObj$cluster_id)

         COPC            EC EpendymalCell     ExDorsal1     ExDorsal2 
           87           332            82           366           354 
    ExDorsal3  ExInVentral1  ExInVentral2   GMAstrocyte      InDorsal 
          149          2382           953          1941           206 
    InVentral    Lymphocyte       Meninge          Mic1          Mic2 
           60           401           208          4479          1900 
           MN         Olig1         Olig2         Olig3         Olig4 
         2481          2228          4294          5847          5126 
        Olig5         Olig6         Olig7         Olig8         Olig9 
         5189          6567          1736           906           336 
         OPC1          OPC2          OPC3      Pericyte   WMAstrocyte 
         2219          1121           110           121          3254 
> table(seuratObj$group_id)

  als   con 
32284 23151 

I really appreciate your help on this! Paria

pariaaliour commented 2 months ago

I also get error when I follow step by step approach:

# Only calculate DE for als condition, with genes that are in the ligand-receptor network
  DE_table <- FindAllMarkers(subset(seuratObj, subset = group_id == "als"),
                             min.pct = 0, logfc.threshold = 0, return.thresh = 1,
                             features = unique(unlist(lr_network_filtered))) 
  Error in h(simpleError(msg, call)) : 
  error in evaluating the argument 'x' in selecting a method for function 'rowSums': subscript out of bounds

Thank you for your help in advance! Paria

csangara commented 2 months ago

Hi,

Can you check that all of the genes in unique(unlist(lr_network_filtered))) exist in your expression data?

Below I show an example of giving random feature names using the example data, which returns a similar error:

> DE_table <- FindAllMarkers(subset(seuratObj, subset = aggregate == "SS"),
                            min.pct = 0, logfc.threshold = 0, return.thresh = 1,
                            features = c("aa", "bb", "cc"))

Calculating cluster CD8 T
Calculating cluster CD4 T
Calculating cluster Treg
Calculating cluster B
Calculating cluster NK
Calculating cluster Mono
Calculating cluster DC
Warning: No DE genes identified
Warning: The following tests were not performed: 
Warning: When testing CD8 T versus all:
    error in evaluating the argument 'x' in selecting a method for function 'rowSums': subscript out of bounds
Warning: When testing CD4 T versus all:
    error in evaluating the argument 'x' in selecting a method for function 'rowSums': subscript out of bounds
Warning: When testing Treg versus all:
    error in evaluating the argument 'x' in selecting a method for function 'rowSums': subscript out of bounds
Warning: When testing B versus all:
    error in evaluating the argument 'x' in selecting a method for function 'rowSums': subscript out of bounds
Warning: When testing NK versus all:
    error in evaluating the argument 'x' in selecting a method for function 'rowSums': subscript out of bounds
Warning: When testing Mono versus all:
    error in evaluating the argument 'x' in selecting a method for function 'rowSums': subscript out of bounds
Warning: When testing DC versus all:
    error in evaluating the argument 'x' in selecting a method for function 'rowSums': subscript out of bounds
pariaaliour commented 2 months ago

Thanks for your reply! When I checked the unique(unlist(lr_network_filtered))I noticed that the column names are printed with genes name. When I removed them it worked. Thank you for your help!