AnnData datatype changes with new version 2.1.0

.github/workflows/python-package.yml

@@ -16,7 +16,7 @@ jobs:
     strategy:
       fail-fast: false
       matrix:
-        python-version: ["3.10", "3.11"]
+        python-version: ["3.11", "3.12"]
 
     steps:
     - uses: actions/checkout@v4

.gitignore

@@ -1,3 +1,5 @@
 test/test_check_with_R_data.py
 rawdata/snRNAseq_cellsWithSubtype.h5ad
 docs/source/notebooks/rawdata/snRNAseq_cellsWithSubtype.h5ad
+docs/source/notebooks/output/HMI_example/Human_Myocardial_Infarction_analysed.h5ad
+tutorials/output/HMI_example/Human_Myocardial_Infarction_analysed.h5ad

.readthedocs.yml

@@ -4,7 +4,7 @@ version: 2
 build:
   os: ubuntu-22.04
   tools:
-    python: "3.10"
+    python: "3.11"
     # You can also specify other tool versions:
     # nodejs: "19"
     # rust: "1.64"

docs/requirements.txt

@@ -7,3 +7,4 @@ sphinx-rtd-size
 gprofiler-official
 sankeyflow
 nbformat
+anndata

docs/source/index.md

@@ -54,7 +54,7 @@ libgeos-dev
 ```
 
 
-## Features v2.0.0
+## Features v2.1.0
 
 
 - Single and Comparative Reports
@@ -64,6 +64,7 @@ libgeos-dev
       - All measures and PC table
       - PC1 and PC2 based barplot
    - LR pair visualization plot can be done
+   - Store analysis results directly in AnnData and export to `.h5ad`
 
 
 # References

pycrosstalker/plots/plot.py

@@ -11,6 +11,7 @@
 from adjustText import adjust_text
 from gprofiler import GProfiler
 from sankeyflow import Sankey
+import json
 
 
 def plot_cci(graph, colors, plt_name, coords, pg, emax=None, leg=False, low=25, high=75, ignore_alpha=False, log=False, efactor=8, vfactor=12, vnames=True, figsize=None, scale_factor=2, node_size=2, font_size=10):
@@ -53,6 +54,11 @@ def plot_cci(graph, colors, plt_name, coords, pg, emax=None, leg=False, low=25,
     Python default plot
 
     """
+    graph = nx.from_pandas_edgelist(graph,
+                                    source='source',
+                                    target='target',
+                                    edge_attr=True,
+                                    create_using=nx.DiGraph())
 
     # Check Maximal Weight
     if emax is None:
@@ -274,14 +280,14 @@ def plot_pca_LR_comparative(lrobj_tblPCA, pca_table, dims=(1, 2), ret=False, ggi
         return pca_plot
 
 
-def plot_bar_rankings(data, table_name, ranking, type = None, filter_sign = None, mode = "cci", top_num = 10):
+def plot_bar_rankings(annData, table_name, ranking, type = None, filter_sign = None, mode = "cci", top_num = 10):
     """
     This function generates the barplot for a given network ranking on the CGI level. Further, the genes can be filtered by selected gene types to filter the plot.
 
     Parameters
     ----------
-    data_object :
-        LRobject with all data
+    annData :
+        AnnData object with all data
 
     table_name :
         name of the ranking table
@@ -302,7 +308,7 @@ def plot_bar_rankings(data, table_name, ranking, type = None, filter_sign = None
     """
 
     if '_x_' in table_name:
-        rankings_table = data['rankings'][table_name]
+        rankings_table = annData.uns['pycrosstalker']['results']['rankings'][table_name]
 
         if type is not None:
             if len(type) == 1:
@@ -820,6 +826,11 @@ def plot_graph_clustermap(graph, weight="LRScore", title="Ligand-Receptor Heatma
     Python Cluster map
 
     """
+    graph = nx.from_pandas_edgelist(graph,
+                                    source='source',
+                                    target='target',
+                                    edge_attr=True,
+                                    create_using=nx.DiGraph())
 
     nodes = list(graph.nodes)
     adj_matrix = nx.to_pandas_adjacency(graph, nodelist=nodes, weight=weight).fillna(0).astype(float)

pycrosstalker/tools/Comparative_condition.py

@@ -1,7 +1,7 @@
 import pandas as pd
 import networkx as nx
 
-def create_diff_table(data, out_path, comparison=None):
+def create_diff_table(annData, out_path, comparison=None):
     """
     Read the lrobject and generate the comparative tables
 
@@ -56,8 +56,11 @@ def process_pair(exp_table, ctr_table):
         freq = (raw_inter - raw_inter.min()) / (raw_inter.max() - raw_inter.min()) + 0.1
         final['freq'] = freq.values
         final['pair'] = final['cellpair']
+        final_data = final_data.drop(final_data.filter(regex='(_x|_y)$').columns, axis=1)
 
         return final_data, final
+
+    data = annData.uns['pycrosstalker']['results']
 
     if comparison is not None:
         for pair in comparison:
@@ -71,13 +74,13 @@ def process_pair(exp_table, ctr_table):
             for _, row in final.iterrows():
                 G.add_edge(row['u'], row['v'], LRScore=row['LRScore'], freq=row['freq'], weight=row['LRScore'], inter=row['freq'])
 
-            data['graphs'][cmp_name] = G
+            data['graphs'][cmp_name] = nx.to_pandas_edgelist(G)
 
             G_ggi = nx.DiGraph()
             for _, row in final_data.iterrows():
                 G_ggi.add_edge(row['ligpair'], row['recpair'], LRScore=row['LRScore'], weight=row['LRScore'], inter=row['LRScore'])
 
-            data['graphs_ggi'][cmp_name] = G_ggi
+            data['graphs_ggi'][cmp_name] = nx.to_pandas_edgelist(G_ggi)
     else:
         ctr_name = list(data['tables'].keys())[0]
         ctr_table = data['tables'][ctr_name]
@@ -91,12 +94,13 @@ def process_pair(exp_table, ctr_table):
             for _, row in final.iterrows():
                 G.add_edge(row['u'], row['v'], LRScore=row['LRScore'], freq=row['freq'], weight=row['LRScore'], inter=row['freq'])
 
-            data['graphs'][cmp_name] = G
+            data['graphs'][cmp_name] = nx.to_pandas_edgelist(G)
 
             G_ggi = nx.DiGraph()
             for _, row in final_data.iterrows():
                 G_ggi.add_edge(row['ligpair'], row['recpair'], LRScore=row['LRScore'], weight=row['LRScore'], inter=row['LRScore'])
 
-            data['graphs_ggi'][cmp_name] = G_ggi
+            data['graphs_ggi'][cmp_name] = nx.to_pandas_edgelist(G_ggi)
 
-    return data
+    annData.uns['pycrosstalker']['results'] = data
+    return annData

pycrosstalker/tools/Single_Condition.py

@@ -4,16 +4,17 @@
 from sklearn.decomposition import PCA
 import matplotlib.pyplot as plt
 import matplotlib.colors as mcolors
+from anndata import AnnData
 from .utils import *
 
-def read_lr_single_condition(lrpaths, sel_columns, out_path="/tmp/", sep=",", colors=None):
+def read_lr_single_condition(input, sel_columns, out_path="/tmp/", sep=",", colors=None):
     """
     This function loads the single conditions LR outputs and use it to generate the report data and it`s object It assumes that the table presents the following columns Ligand, Ligand.Cluster, Receptor,Receptor.Cluster and MeanLR/another measure
 
     Parameters
     ----------
-    lrpaths :
-        Named vector with the lrpaths of each output
+    input :
+        Named vector with the lrpaths of each output or an AnnData object
     sel_columns :
         selected columns
     out_path :
@@ -28,21 +29,27 @@ def read_lr_single_condition(lrpaths, sel_columns, out_path="/tmp/", sep=",", co
     LRObject
 
     """
+    if isinstance(input, dict):
+        annData = AnnData(uns = {'pycrosstalker': {'path': {}}})
+        for cond, lrpath in input.items():
+            # Reading data
+            if isinstance(lrpath, str):
+                data1 = pd.read_csv(lrpath, sep=sep)
+                annData.uns['pycrosstalker']['path'][cond] = data1
+            else:
+                raise ValueError("Issue with input paths, please check!")
+    elif isinstance(input, AnnData):
+        annData = input
+    else:
+        raise ValueError("Input parameter must be either a file path or an AnnData object")
 
     data = {}
     graphs = {}
     graphs_ggi = {}
     unif_celltypes = []
 
-    for cond, lrpath in lrpaths.items():
-        # Reading data
-        if isinstance(lrpath, str):
-            data1 = pd.read_csv(lrpath, sep=sep)    
-        elif isinstance(lrpath, pd.DataFrame):
-            data1 = lrpath.copy()
-        else:
-            raise ValueError("lrpath must be either a file path or a DataFrame")
-
+    for cond in annData.uns['pycrosstalker']['path'].keys():
+        data1 = annData.uns['pycrosstalker']['path'][cond]
         if not (data1['gene_A'].str.contains(r'\|').sum() > 0):
             data1 = add_node_type(data1)
 
@@ -97,8 +104,8 @@ def read_lr_single_condition(lrpaths, sel_columns, out_path="/tmp/", sep=",", co
                             inter=row['LRScore'])  # Add thickness (inter)
 
         data[cond] = data1
-        graphs[cond] = graph1
-        graphs_ggi[cond] = graph2
+        graphs[cond] = nx.to_pandas_edgelist(graph1)
+        graphs_ggi[cond] = nx.to_pandas_edgelist(graph2)
 
     # Create a full graph
     template = nx.complete_graph(len(set(unif_celltypes)))
@@ -119,4 +126,5 @@ def read_lr_single_condition(lrpaths, sel_columns, out_path="/tmp/", sep=",", co
           "pca" : {},
           "stats" : {}} 
 
-    return lr
+    annData.uns['pycrosstalker']['results'] = lr
+    return annData

pycrosstalker/tools/generate_report.py

@@ -1,12 +1,13 @@
 import os
 import pandas as pd
 from tqdm import tqdm
+from anndata import AnnData
 from .Comparative_condition import *
 from .Single_Condition import *
 from .utils import *
 
 
-def analise_LR(lrpaths, 
+def analise_LR(input,
                 genes=None, 
                 tf_genes=None, 
                 out_path=None, 
@@ -16,15 +17,15 @@ def analise_LR(lrpaths,
                 out_file=None, 
                 output_fmt="html_document", 
                 sel_columns=['source','target','gene_A','gene_B','type_gene_A','type_gene_B','MeanLR'], 
-                org='hsa', comparison=None, filtered_net=False):
+                org='hsa', comparison=None, filtered_net=False, filename=None):
 
     """
     Core engine to generate report. Here we perform all the computation related to pyCrossTalkeR
 
     Parameters
     ----------
-    lrpaths :
-        Paths of single condition LR data
+    input :
+        Named vector with the lrpaths of each output or an AnnData object    
     genes :
         list of genes to be considered in the sankey plots
     out_path :
@@ -41,34 +42,40 @@ def analise_LR(lrpaths,
         rmarkdown render output format parameter
     sel_columns :
         columns from data
+    filename :
+        filename prefix for output files (to be provided if not already present in AnnData object)
 
     Returns
     -------
     Rmarkdown report all objects from each step
 
     """
 
-    data = read_lr_single_condition(lrpaths, 
+    annData = read_lr_single_condition(input,
                                     sel_columns, 
                                     out_path, 
                                     sep, 
                                     colors)
 
+    if isinstance(input, dict) or (isinstance(input, AnnData) and 'details' not in annData.uns['pycrosstalker']):
+        annData.uns['pycrosstalker']['details'] = {}
+        annData.uns['pycrosstalker']['details']['filename'] = filename
+
     print("Create a Differential Table")
-    if len(lrpaths) > 1:
-        data = create_diff_table(data, out_path, comparison)
-        data = fisher_test_cci(data, 'LRScore', out_path, comparison)
-        data = mannwhitneyu_test_cci(data, 'LRScore', out_path, comparison)
-        data = filtered_graphs(data, out_path)
+    if len(annData.uns['pycrosstalker']['path']) > 1:
+        annData = create_diff_table(annData, out_path, comparison)
+        annData = fisher_test_cci(annData, 'LRScore', out_path, comparison)
+        annData = mannwhitneyu_test_cci(annData, 'LRScore', out_path, comparison)
+        annData = filtered_graphs(annData, out_path)
 
     print("Calculating CCI Ranking")
-    data = ranking(data, out_path, sel_columns=sel_columns, slot="graphs")
+    annData = ranking(annData, out_path, sel_columns=sel_columns, slot="graphs")
     print("Calculating GCI Ranking")
-    data = ranking(data, out_path, sel_columns=sel_columns, slot="graphs_ggi")
+    annData = ranking(annData, out_path, sel_columns=sel_columns, slot="graphs_ggi")
     print("Network Analysis Done")
 
-    with open(os.path.join(out_path, "LR_data.pkl"), "wb") as f:
-        pickle.dump(data, f)
+    print("Generating h5ad file with Analysed Results")
+    annData.write(os.path.join(out_path, annData.uns['pycrosstalker']['details']['filename'] + "_analysed.h5ad"))
 
-    return(data)
+    return(annData)