search

KellisLab / scdemon

single-cell decorrelated module networks
https://scdemon.readthedocs.io/
GNU General Public License v3.0
18 stars 0 forks source link

Tutorial does not work. #2

Open jhaberbe opened 2 months ago

jhaberbe commented 2 months ago

Hi,

I was interested in getting your tutorial to work, but ran into some issues. When I ran the code on you github repo:

import numpy as np
import scanpy as sc
import scdemon as sm
from scdemon.utils import recipe_full
from scdemon import plotting as pl

# Load one of scanpy's datasets:
adata = sc.datasets.pbmc3k()
recipe_full(adata, preprocess=True, annotate=True)

# Make the modules handling object:
mod = sm.modules_core(adata, suffix='pbmc_example', k=100)
mod.setup()

# Build a gene-graph graph and learn modules:
graph_id = 'base'
mod.make_graph(graph_id, resolution=2.5)

# Plot genes on the gene-gene graph and on the gene-level UMAP basis
pl.plot_genes(mod, graph_id, attr="leiden", show_labels=True, width=16)
pl.plot_genes(mod, graph_id, basis='umap', attr="leiden", width=16)

# Plot module expression on the cell-level UMAP basis:
pl.plot_umap_grid(mod, graph_id)

# Get the modules and/or print them out:
mlist = mod.get_modules(graph_id, print_modules=False)
mod.save_modules(graph_id)

# Get functional enrichments for the modules:
gpres = sm.get_goterms(mod, graph_id)

I recieved the following error

---------------------------------------------------------------------------
AttributeError                            Traceback (most recent call last)
Cell In[51], line 12
      9 recipe_full(adata, preprocess=True, annotate=True)
     11 # Make the modules handling object:
---> 12 mod = sm.modules_core(adata, suffix='pbmc_example', k=100)
     13 mod.setup()
     15 # Build a gene-graph graph and learn modules:

AttributeError: module 'scdemon' has no attribute 'modules_core'

While this works after I change the scdemon argument from 'sm.modules_core' to 'sm.modules' for your test scenario, I run into the following problem attempting to use my own dataset:

import numpy as np
import scanpy as sc
import scdemon as sm
from scdemon.utils import recipe_full
from scdemon import plotting as pl

# Load one of scanpy's datasets:
adata = sc.read_h5ad("path/to/my/data.h5ad")
mod = sm.modules(adata, suffix='pbmc_example', k=100)

mod.setup()

# Build a gene-graph graph and learn modules:
graph_id = 'base'
mod.make_graph(graph_id, resolution=2.5)

# Plot genes on the gene-gene graph and on the gene-level UMAP basis
pl.plot_genes(mod, graph_id, attr="leiden", show_labels=True, width=16)
pl.plot_genes(mod, graph_id, basis='umap', attr="leiden", width=16)

# Plot module expression on the cell-level UMAP basis:
pl.plot_umap_grid(mod, graph_id)

# Get the modules and/or print them out:
mlist = mod.get_modules(graph_id, print_modules=False)
mod.save_modules(graph_id)

# Get functional enrichments for the modules:
gpres = sm.get_goterms(mod, graph_id)
---------------------------------------------------------------------------
ValueError                                Traceback (most recent call last)
Cell In[58], line 5
      3 # Build a gene-graph graph and learn modules:
      4 graph_id = 'base'
----> 5 mod.make_graph(graph_id, resolution=2.5)
      7 # Plot genes on the gene-gene graph and on the gene-level UMAP basis
      8 pl.plot_genes(mod, graph_id, attr="leiden", show_labels=True, width=16)

File [/oak/stanford/projects/kibr/Reorganizing/Projects/choroid-plexus/.venv/lib/python3.12/site-packages/scdemon/modules.py:283](https://vscode-remote+ondemand-002esherlock-002estanford-002eedu.vscode-resource.vscode-cdn.net/oak/stanford/projects/kibr/Reorganizing/Projects/choroid-plexus/.venv/lib/python3.12/site-packages/scdemon/modules.py:283), in modules.make_graph(self, graph_id, multigraph, power, **kwargs)
    281     self._make_merged_graph(graph_id, power_list=power, **kwargs)
    282 else:
--> 283     self._make_single_graph(graph_id, power=power, **kwargs)

File [/oak/stanford/projects/kibr/Reorganizing/Projects/choroid-plexus/.venv/lib/python3.12/site-packages/scdemon/modules.py:322](https://vscode-remote+ondemand-002esherlock-002estanford-002eedu.vscode-resource.vscode-cdn.net/oak/stanford/projects/kibr/Reorganizing/Projects/choroid-plexus/.venv/lib/python3.12/site-packages/scdemon/modules.py:322), in modules._make_single_graph(self, graph_id, power, filter_covariate, raw, method, resolution, adjacency_only, full_graph_only, keep_all_z, layout, **kwargs)
    318     self.graphs[graph_id].construct_graph(
    319         full_graph=True, modules=False, layout=False)
    320 else:
    321     # Build the graph, get modules, and annotate genes:
--> 322     self.graphs[graph_id].construct_graph(
    323         resolution=resolution, layout=layout)
    324     self.graphs[graph_id].populate_modules(self.adata, attr='leiden')
    325     self.graphs[graph_id].match_genes_to_modules(attr='leiden')

File [/oak/stanford/projects/kibr/Reorganizing/Projects/choroid-plexus/.venv/lib/python3.12/site-packages/scdemon/graph/gene_graph.py:105](https://vscode-remote+ondemand-002esherlock-002estanford-002eedu.vscode-resource.vscode-cdn.net/oak/stanford/projects/kibr/Reorganizing/Projects/choroid-plexus/.venv/lib/python3.12/site-packages/scdemon/graph/gene_graph.py:105), in gene_graph.construct_graph(self, resolution, method, full_graph, modules, layout)
    102 # NOTE: check if modules is None, also remake if made new graph
    103 if modules:
    104     # 7. Perform community detection (NMF, BigClam, Leiden)
--> 105     self.calculate_gene_modules(
    106         resolution=resolution, method=method)
    107 if layout:
    108     # 7a. Layout graph
    109     self.layout_graph(layout_method=self.layout_method)

File [/oak/stanford/projects/kibr/Reorganizing/Projects/choroid-plexus/.venv/lib/python3.12/site-packages/scdemon/graph/gene_graph.py:185](https://vscode-remote+ondemand-002esherlock-002estanford-002eedu.vscode-resource.vscode-cdn.net/oak/stanford/projects/kibr/Reorganizing/Projects/choroid-plexus/.venv/lib/python3.12/site-packages/scdemon/graph/gene_graph.py:185), in gene_graph.calculate_gene_modules(self, method, **kwargs)
    181 else:
    182     # NOTE: Could also implement louvain, HDBSCAN, etc
    183     raise ValueError("Method '%s' not found" % method)
--> 185 self.colors[method] = set_module_colors(
    186     self.assign[method], self.snapcols)

File [/oak/stanford/projects/kibr/Reorganizing/Projects/choroid-plexus/.venv/lib/python3.12/site-packages/scdemon/graph/utils_community_detection.py:68](https://vscode-remote+ondemand-002esherlock-002estanford-002eedu.vscode-resource.vscode-cdn.net/oak/stanford/projects/kibr/Reorganizing/Projects/choroid-plexus/.venv/lib/python3.12/site-packages/scdemon/graph/utils_community_detection.py:68), in set_module_colors(assign, palette)
     65 def set_module_colors(assign, palette):
     66     # Assign colors to each cluster:
     67     # TODO: Fix so it works for overlapping modules too
---> 68     nclust = np.max(assign) + 1
     69     ncols = len(palette)
     70     rep_col = int(np.ceil((nclust + 2) [/](https://vscode-remote+ondemand-002esherlock-002estanford-002eedu.vscode-resource.vscode-cdn.net/) (ncols * 1.0)))

File [/oak/stanford/projects/kibr/Reorganizing/Projects/choroid-plexus/.venv/lib/python3.12/site-packages/numpy/core/fromnumeric.py:2810](https://vscode-remote+ondemand-002esherlock-002estanford-002eedu.vscode-resource.vscode-cdn.net/oak/stanford/projects/kibr/Reorganizing/Projects/choroid-plexus/.venv/lib/python3.12/site-packages/numpy/core/fromnumeric.py:2810), in max(a, axis, out, keepdims, initial, where)
   2692 @array_function_dispatch(_max_dispatcher)
   2693 @set_module('numpy')
   2694 def max(a, axis=None, out=None, keepdims=np._NoValue, initial=np._NoValue,
   2695          where=np._NoValue):
   2696     """
   2697     Return the maximum of an array or maximum along an axis.
   2698 
   (...)
   2808     5
   2809     """
-> 2810     return _wrapreduction(a, np.maximum, 'max', axis, None, out,
   2811                           keepdims=keepdims, initial=initial, where=where)

File [/oak/stanford/projects/kibr/Reorganizing/Projects/choroid-plexus/.venv/lib/python3.12/site-packages/numpy/core/fromnumeric.py:88](https://vscode-remote+ondemand-002esherlock-002estanford-002eedu.vscode-resource.vscode-cdn.net/oak/stanford/projects/kibr/Reorganizing/Projects/choroid-plexus/.venv/lib/python3.12/site-packages/numpy/core/fromnumeric.py:88), in _wrapreduction(obj, ufunc, method, axis, dtype, out, **kwargs)
     85         else:
     86             return reduction(axis=axis, out=out, **passkwargs)
---> 88 return ufunc.reduce(obj, axis, dtype, out, **passkwargs)

ValueError: zero-size array to reduction operation maximum which has no identity

I'm using python/3.12.1

cboix commented 2 months ago

Hi,

Thanks for catching this. I fixed the issue with the tutorial, both in the readme and in the documentation.

I am not sure why your dataset is not working. It seems like the code fails when trying to assign modules to genes. It may be the case that no genes have gotten modules assigned or that there is only one module.

Can you check what mod.assign looks like? mod.assign['leiden'] should be a list of the assignments of modules to genes. From the error, it looks like it's empty.

One possibility is that too many genes were filtered out or that the z-score cutoff for calling gene-gene edges is too high. One thing that you can try is to lower the cutoff, by passing z=3 or lower (default is z=4.5) to make_graph. Let me know if that works, I can update documentation to match. I know it isn't super clear at the moment, we're still working on refactoring the code to make it cleaner.