Repository: aertslab/pySCENIC (language: python)

Issue:
pyscenic ctx - ValueError: The columns in the computed data do not match the columns in the provided metadata
**Describe the bug**
Hi, I got a ValueError when runninig "pyscenic ctx". I follow the Nature Protocol and download the ranking database in feather format. The list of TF I used is SCENICprotocol/example/allTFs_hg38.txt, and the gene name in my matrix is also Gene Symbol. I found that some TFs in the adj.csv are not in the feather files (colnames), but after I excluded them, I still got the error. Could you please help me? 

**Steps to reproduce the behavior**
1. Command run when the error occurred:
<!-- Please specify the command used (if applicable, otherwise delete this block): -->
```bash
pyscenic grn $f_loom_path SCENICprotocol/example/allTFs_hg38.txt-o adj.csv --num_workers 20 --seed 777

pyscenic ctx ./scenic_results/adj.csv /BGFS1/home/aic/tools/scenic/resourc
es/hg38__refseq-r80__10kb_up_and_down_tss.mc9nr.feather \
> --annotations_fname ~/tools/scenic/resources/motifs-v9-nr.hgnc-m0.001-o0.0.tbl \
> --expression_mtx_fname /BGFS1/home/aic/projects/FTE/submission04/results/20200721_cancerA/CD4Treg_cancerA_sub13467_clustered_annotated_C14.loom \
> --mode "dask_multiprocessing" \
> --output ./scenic_results/reg.csv \
> --num_workers 20 \
>  --mask_dropouts

2020-08-11 10:13:20,532 - pyscenic.cli.pyscenic - INFO - Creating modules.

2020-08-11 10:13:22,680 - pyscenic.cli.pyscenic - INFO - Loading expression matrix.

2020-08-11 10:13:29,778 - pyscenic.utils - INFO - Calculating Pearson correlations.

2020-08-11 10:13:29,778 - pyscenic.utils - WARNING - Note on correlation calculation: the default behaviour for calculating the correlations has changed after pySCENIC verion 0.9.16. Previously, the default was to calculate the correlation between a TF and target gene using only cells with non-zero expression values (mask_dropouts=True). The current default is now to use all cells to match the behavior of the R verision of SCENIC. The original settings can be retained by setting 'rho_mask_dropouts=True' in the modules_from_adjacencies function, or '--mask_dropouts' from the CLI.
        Dropout masking is currently set to [True].

/BGFS1/home/aic/.conda/envs/sc/lib/python3.8/site-packages/pyscenic/utils.py:138: RuntimeWarning: invalid value encountered in greater
  regulations = (rhos > rho_threshold).astype(int) - (rhos < -rho_threshold).astype(int)
/BGFS1/home/aic/.conda/envs/sc/lib/python3.8/site-packages/pyscenic/utils.py:138: RuntimeWarning: invalid value encountered in less
  regulations = (rhos > rho_threshold).astype(int) - (rhos < -rho_threshold).astype(int)

2020-08-11 10:19:55,129 - pyscenic.utils - INFO - Creating modules.

2020-08-11 10:20:50,275 - pyscenic.cli.pyscenic - INFO - Loading databases.

2020-08-11 10:20:50,276 - pyscenic.cli.pyscenic - INFO - Calculating regulons.
[                                        ] | 0% Completed | 46.6s
Traceback (most recent call last):
  File "/BGFS1/home/aic/.conda/envs/sc/bin/pyscenic", line 8, in <module>
    sys.exit(main())
  File "/BGFS1/home/aic/.conda/envs/sc/lib/python3.8/site-packages/pyscenic/cli/pyscenic.py", line 420, in main
    args.func(args)
  File "/BGFS1/home/aic/.conda/envs/sc/lib/python3.8/site-packages/pyscenic/cli/pyscenic.py", line 159, in prune_targets_command
    df_motifs = calc_func(dbs, modules, motif_annotations_fname,
  File "/BGFS1/home/aic/.conda/envs/sc/lib/python3.8/site-packages/pyscenic/prune.py", line 349, in prune2df
    return _distributed_calc(rnkdbs, modules, motif_annotations_fname, transformation_func, aggregation_func,
  File "/BGFS1/home/aic/.conda/envs/sc/lib/python3.8/site-packages/pyscenic/prune.py", line 300, in _distributed_calc
    return create_graph().compute(scheduler='processes', num_workers=num_workers if num_workers else cpu_count())
  File "/BGFS1/home/aic/.conda/envs/sc/lib/python3.8/site-packages/dask/base.py", line 166, in compute
    (result,) = compute(self, traverse=False, **kwargs)
  File "/BGFS1/home/aic/.conda/envs/sc/lib/python3.8/site-packages/dask/base.py", line 437, in compute
    results = schedule(dsk, keys, **kwargs)
  File "/BGFS1/home/aic/.conda/envs/sc/lib/python3.8/site-packages/dask/multiprocessing.py", line 212, in get
    result = get_async(
  File "/BGFS1/home/aic/.conda/envs/sc/lib/python3.8/site-packages/dask/local.py", line 486, in get_async
    raise_exception(exc, tb)
  File "/BGFS1/home/aic/.conda/envs/sc/lib/python3.8/site-packages/dask/local.py", line 316, in reraise
    raise exc
  File "/BGFS1/home/aic/.conda/envs/sc/lib/python3.8/site-packages/dask/local.py", line 222, in execute_task
    result = _execute_task(task, data)
  File "/BGFS1/home/aic/.conda/envs/sc/lib/python3.8/site-packages/dask/core.py", line 119, in _execute_task
    return func(*args2)
  File "/BGFS1/home/aic/.conda/envs/sc/lib/python3.8/site-packages/dask/dataframe/utils.py", line 653, in check_meta
    check_matching_columns(meta, x)
  File "/BGFS1/home/aic/.conda/envs/sc/lib/python3.8/site-packages/dask/dataframe/utils.py", line 674, in check_matching_columns
    raise ValueError(
ValueError: The columns in the computed data do not match the columns in the provided metadata
  Extra:   []
  Missing: []
```

**Please complete the following information:**
- pySCENIC version: 0.10.0
- Installation method: Pip
- Run environment:  HPC
- OS: (Ubuntu 7.4.0-1ubuntu1~18.04.1))
- Package versions: [obtain using `pip freeze`, `conda list`, or skip this if using Docker/Singularity]:
<!-- Put your package version list in this code block (if applicable, else delete the block): -->
```
# Name                    Version                   Build  Channel
_libgcc_mutex             0.1                 conda_forge    conda-forge
_openmp_mutex             4.5                      1_llvm    conda-forge
anndata                   0.7.1                    pypi_0    pypi
arboreto                  0.1.5                    pypi_0    pypi
attrs                     19.3.0                     py_0    conda-forge
backcall                  0.1.0                      py_0    conda-forge
blas                      1.0                    openblas
bleach                    3.1.5              pyh9f0ad1d_0    conda-forge
blosc                     1.18.1               he1b5a44_0    conda-forge
boltons                   20.2.0                   pypi_0    pypi
bzip2                     1.0.8                h516909a_2    conda-forge
ca-certificates           2020.6.24                     0
cairo                     1.16.0            hcf35c78_1003    conda-forge
certifi                   2020.6.20                py38_0
cffi                      1.14.0                   pypi_0    pypi
chardet                   3.0.4                    pypi_0    pypi
click                     7.1.2                    pypi_0    pypi
cloudpickle               1.5.0                    pypi_0    pypi
cmake                     3.17.0               h28c56e5_0    conda-forge
cycler                    0.10.0                     py_2    conda-forge
cytoolz                   0.10.1                   pypi_0    pypi
dask                      2.11.0                   pypi_0    pypi
dbus                      1.13.6               he372182_0    conda-forge
decorator                 4.4.2                      py_0    conda-forge
defusedxml                0.6.0                      py_0    conda-forge
dill                      0.3.2                    pypi_0    pypi
distributed               2.11.0                   pypi_0    pypi
entrypoints               0.3             py38h32f6830_1001    conda-forge
expat                     2.2.9                he1b5a44_2    conda-forge
fa2                       0.3.5                    pypi_0    pypi
fontconfig                2.13.1            h86ecdb6_1001    conda-forge
freetype                  2.10.1               he06d7ca_0    conda-forge
frozendict                1.2                      pypi_0    pypi
fsspec                    0.8.0                    pypi_0    pypi
get-version               2.1                      pypi_0    pypi
gettext                   0.19.8.1          hc5be6a0_1002    conda-forge
glib                      2.64.2               h6f030ca_0    conda-forge
gmp                       6.2.0                he1b5a44_2    conda-forge
gprofiler-official        1.0.0                    pypi_0    pypi
gsl                       2.6                  h294904e_0    conda-forge
gst-plugins-base          1.14.5               h0935bb2_2    conda-forge
gstreamer                 1.14.5               h36ae1b5_2    conda-forge
h5py                      2.10.0          nompi_py38h513d04c_102    conda-forge
hdf5                      1.10.5          nompi_h3c11f04_1104    conda-forge
heapdict                  1.0.1                    pypi_0    pypi
icu                       64.2                 he1b5a44_1    conda-forge
idna                      2.9                      pypi_0    pypi
igraph                    0.7.1             h9e3b1fc_1007    conda-forge
importlib-metadata        1.6.0            py38h32f6830_0    conda-forge
importlib_metadata        1.6.0                         0    conda-forge
interlap                  0.2.6                    pypi_0    pypi
ipykernel                 5.2.1            py38h23f93f0_0    conda-forge
ipython                   7.14.0           py38h32f6830_0    conda-forge
ipython_genutils          0.2.0                      py_1    conda-forge
ipywidgets                7.5.1                      py_0    conda-forge
jedi                      0.17.0           py38h32f6830_0    conda-forge
jinja2                    2.11.2             pyh9f0ad1d_0    conda-forge
joblib                    0.14.1                     py_0    conda-forge
jpeg                      9c                h14c3975_1001    conda-forge
jsonschema                3.2.0            py38h32f6830_1    conda-forge
jupyter                   1.0.0                      py_2    conda-forge
jupyter_client            6.1.3                      py_0    conda-forge
jupyter_console           6.1.0                      py_1    conda-forge
jupyter_core              4.6.3            py38h32f6830_1    conda-forge
kiwisolver                1.2.0            py38hbf85e49_0    conda-forge
krb5                      1.17.1               h2fd8d38_0    conda-forge
ld_impl_linux-64          2.34                 h53a641e_0    conda-forge
legacy-api-wrap           1.2                      pypi_0    pypi
libblas                   3.8.0               16_openblas    conda-forge
libcblas                  3.8.0               16_openblas    conda-forge
libclang                  9.0.1           default_hde54327_0    conda-forge
libcurl                   7.69.1               hf7181ac_0    conda-forge
libedit                   3.1.20170329      hf8c457e_1001    conda-forge
libffi                    3.2.1             he1b5a44_1007    conda-forge
libgcc-ng                 9.2.0                h24d8f2e_2    conda-forge
libgfortran-ng            7.3.0                hdf63c60_5    conda-forge
libiconv                  1.15              h516909a_1006    conda-forge
liblapack                 3.8.0               16_openblas    conda-forge
libllvm8                  8.0.1                hc9558a2_0    conda-forge
libllvm9                  9.0.1                he513fc3_1    conda-forge
libopenblas               0.3.9                h5ec1e0e_0    conda-forge
libpng                    1.6.37               hed695b0_1    conda-forge
libsodium                 1.0.17               h516909a_0    conda-forge
libssh2                   1.8.2                h22169c7_2    conda-forge
libstdcxx-ng              9.2.0                hdf63c60_2    conda-forge
libuuid                   2.32.1            h14c3975_1000    conda-forge
libuv                     1.34.0               h516909a_0    conda-forge
libxcb                    1.13              h14c3975_1002    conda-forge
libxkbcommon              0.10.0               he1b5a44_0    conda-forge
libxml2                   2.9.10               hee79883_0    conda-forge
llvm-openmp               10.0.0               hc9558a2_0    conda-forge
llvmlite                  0.31.0           py38h4f45e52_1    conda-forge
loompy                    3.0.6                    pypi_0    pypi
louvain                   0.6.1            py38he1b5a44_2    conda-forge
lzo                       2.10              h14c3975_1000    conda-forge
markupsafe                1.1.1            py38h1e0a361_1    conda-forge
matplotlib-base           3.2.1            py38h2af1d28_0    conda-forge
mistune                   0.8.4           py38h1e0a361_1001    conda-forge
mock                      4.0.2            py38h32f6830_0    conda-forge
msgpack                   1.0.0                    pypi_0    pypi
multicoretsne             0.1                      pypi_0    pypi
multiprocessing-on-dill   3.5.0a4                  pypi_0    pypi
natsort                   7.0.1                      py_0    conda-forge
nbconvert                 5.6.1            py38h32f6830_1    conda-forge
nbformat                  5.0.6                      py_0    conda-forge
ncurses                   6.1               hf484d3e_1002    conda-forge
networkx                  2.4                        py_1    conda-forge
notebook                  6.0.3                    py38_0    conda-forge
nspr                      4.25                 he1b5a44_0    conda-forge
nss                       3.47                 he751ad9_0    conda-forge
numba                     0.48.0           py38hb3f55d8_0    conda-forge
numexpr                   2.7.1            py38hcb8c335_1    conda-forge
numpy                     1.18.4           py38h8854b6b_0    conda-forge
numpy-groupies            0+unknown                pypi_0    pypi
openssl                   1.1.1g               h7b6447c_0
packaging                 20.1                       py_0    conda-forge
pandas                    0.25.3                   pypi_0    pypi
pandoc                    2.9.2.1                       0    conda-forge
pandocfilters             1.4.2                      py_1    conda-forge
parso                     0.7.0              pyh9f0ad1d_0    conda-forge
patsy                     0.5.1                      py_0    conda-forge
pcre                      8.44                 he1b5a44_0    conda-forge
pexpect                   4.8.0            py38h32f6830_1    conda-forge
pickleshare               0.7.5           py38h32f6830_1001    conda-forge
pip                       20.1               pyh9f0ad1d_0    conda-forge
pixman                    0.38.0            h516909a_1003    conda-forge
prometheus_client         0.7.1                      py_0    conda-forge
prompt-toolkit            3.0.5                      py_0    conda-forge
prompt_toolkit            3.0.5                         0    conda-forge
psutil                    5.7.2                    pypi_0    pypi
pthread-stubs             0.4               h14c3975_1001    conda-forge
ptyprocess                0.6.0                   py_1001    conda-forge
pyarrow                   0.16.0                   pypi_0    pypi
pycairo                   1.19.1           py38h323dad1_3    conda-forge
pycparser                 2.20                     pypi_0    pypi
pygments                  2.6.1                      py_0    conda-forge
pyparsing                 2.4.7              pyh9f0ad1d_0    conda-forge
pyqt                      5.12.3           py38ha8c2ead_3    conda-forge
pyqt5-sip                 4.19.18                  pypi_0    pypi
pyqtchart                 5.12                     pypi_0    pypi
pyqtwebengine             5.12.1                   pypi_0    pypi
pyrsistent                0.16.0           py38h1e0a361_0    conda-forge
pyscenic                  0.10.0                   pypi_0    pypi
pytables                  3.6.1            py38h9f153d1_1    conda-forge
python                    3.8.2           he5300dc_7_cpython    conda-forge
python-dateutil           2.8.1                      py_0    conda-forge
python-igraph             0.8.2            py38hb721cf7_0    conda-forge
python_abi                3.8                      1_cp38    conda-forge
pytz                      2020.1             pyh9f0ad1d_0    conda-forge
pyyaml                    5.3.1                    pypi_0    pypi
pyzmq                     19.0.0           py38ha71036d_1    conda-forge
qt                        5.12.5               hd8c4c69_1    conda-forge
qtconsole                 4.7.3              pyh9f0ad1d_0    conda-forge
qtpy                      1.9.0                      py_0    conda-forge
readline                  8.0                  hf8c457e_0    conda-forge
requests                  2.23.0                   pypi_0    pypi
rhash                     1.3.6             h14c3975_1001    conda-forge
scanpy                    1.4.6                    pypi_0    pypi
schpf                     0.3.0                    pypi_0    pypi
scikit-learn              0.23.1           py38h7ea95a0_0
scipy                     1.4.1            py38h18bccfc_3    conda-forge
seaborn                   0.10.1                     py_0    conda-forge
send2trash                1.5.0                      py_0    conda-forge
setuptools                46.1.3           py38h32f6830_0    conda-forge
setuptools-scm            3.5.0                    pypi_0    pypi
simplegeneric             0.8.1                      py_1    conda-forge
six                       1.14.0                     py_1    conda-forge
sortedcontainers          2.2.2                    pypi_0    pypi
sqlite                    3.30.1               hcee41ef_0    conda-forge
statsmodels               0.11.1           py38h7b6447c_0
tbb                       2020.0.133               pypi_0    pypi
tblib                     1.7.0                    pypi_0    pypi
terminado                 0.8.3            py38h32f6830_1    conda-forge
testpath                  0.4.4                      py_0    conda-forge
texttable                 1.6.2                      py_0    conda-forge
threadpoolctl             2.1.0              pyh5ca1d4c_0
tk                        8.6.10               hed695b0_0    conda-forge
toolz                     0.10.0                   pypi_0    pypi
tornado                   6.0.4            py38h1e0a361_1    conda-forge
tqdm                      4.46.0             pyh9f0ad1d_0    conda-forge
traitlets                 4.3.3            py38h32f6830_1    conda-forge
tzlocal                   2.0.0                      py_0    conda-forge
umap-learn                0.4.2                    pypi_0    pypi
urllib3                   1.25.9                   pypi_0    pypi
wcwidth                   0.1.9              pyh9f0ad1d_0    conda-forge
webencodings              0.5.1                      py_1    conda-forge
wheel                     0.34.2                     py_1    conda-forge
widgetsnbextension        3.5.1                    py38_0    conda-forge
xorg-kbproto              1.0.7             h14c3975_1002    conda-forge
xorg-libice               1.0.10               h516909a_0    conda-forge
xorg-libsm                1.2.3             h84519dc_1000    conda-forge
xorg-libx11               1.6.9                h516909a_0    conda-forge
xorg-libxau               1.0.9                h14c3975_0    conda-forge
xorg-libxdmcp             1.1.3                h516909a_0    conda-forge
xorg-libxext              1.3.4                h516909a_0    conda-forge
xorg-libxrender           0.9.10            h516909a_1002    conda-forge
xorg-renderproto          0.11.1            h14c3975_1002    conda-forge
xorg-xextproto            7.3.0             h14c3975_1002    conda-forge
xorg-xproto               7.0.31            h14c3975_1007    conda-forge
xz                        5.2.5                h516909a_0    conda-forge
zeromq                    4.3.2                he1b5a44_2    conda-forge
zict                      2.0.0                    pypi_0    pypi
zipp                      3.1.0                      py_0    conda-forge
zlib                      1.2.11            h516909a_1006    conda-forge
```

AttributeError: module 'numba' has no attribute 'jitclass'
Dear friends,
When I performed the pyscenic analysis, it show the error as follows,
In [159]: !pyscenic grn {f_loom_path_scenic} {f_tfs} -o adj.csv --num_workers 20
/home/lee/miniconda3/lib/python3.8/site-packages/dask/config.py:161: YAMLLoadWarning: calling yaml.load() without Loader=... is deprecated, as the default Loader is unsafe. Please read https://msg.pyyaml.org/load for full details.
  data = yaml.load(f.read()) or {}
Traceback (most recent call last):
  File "/home/lee/miniconda3/bin/pyscenic", line 5, in <module>
    from pyscenic.cli.pyscenic import main
  File "/home/lee/miniconda3/lib/python3.8/site-packages/pyscenic/cli/pyscenic.py", line 19, in <module>
    from pyscenic.utils import modules_from_adjacencies, add_correlation
  File "/home/lee/miniconda3/lib/python3.8/site-packages/pyscenic/utils.py", line 6, in <module>
    from .math import masked_rho4pairs
  File "/home/lee/miniconda3/lib/python3.8/site-packages/pyscenic/math.py", line 5, in <module>
    from numba import *
AttributeError: module 'numba' has no attribute 'jitclass'

Do you know how to solve this problem? I try to solve this problem through the idea from https://github.com/Roth-Lab/pyclone-vi/issues/7. But it still doesn't work.
Thanks
Error when running modules = list(modules_from_adjacencies(adjacencies, ex_matrix))
Hello, 

I am getting this error the running modules = list(modules_from_adjacencies(adjacencies, ex_matrix)) : 

2020-03-23 13:46:47,461 - pyscenic.utils - INFO - Calculating Pearson correlations.

2020-03-23 13:46:47,461 - pyscenic.utils - WARNING - Note on correlation calculation: the default behaviour for calculating the correlations has changed after pySCENIC verion 0.9.16. Previously, the default was to calculate the correlation between a TF and target gene using only cells with non-zero expression values (mask_dropouts=True). The current default is now to use all cells to match the behavior of the R verision of SCENIC. The original settings can be retained by setting 'rho_mask_dropouts=True' in the modules_from_adjacencies function, or '--mask_dropouts' from the CLI.
	Dropout masking is currently set to [False].
Traceback (most recent call last):
  File "/home/Morgane/anaconda3/lib/python3.7/site-packages/pandas/core/indexes/base.py", line 2897, in get_loc
    return self._engine.get_loc(key)
  File "pandas/_libs/index.pyx", line 107, in pandas._libs.index.IndexEngine.get_loc
  File "pandas/_libs/index.pyx", line 128, in pandas._libs.index.IndexEngine.get_loc
  File "pandas/_libs/index_class_helper.pxi", line 91, in pandas._libs.index.Int64Engine._check_type
KeyError: 'RPS19'

During handling of the above exception, another exception occurred:

Traceback (most recent call last):
  File "<stdin>", line 1, in <module>
  File "/home/Morgane/anaconda3/lib/python3.7/site-packages/pyscenic/utils.py", line 268, in modules_from_adjacencies
    rho_threshold=rho_threshold, mask_dropouts=rho_mask_dropouts)
  File "/home/Morgane/anaconda3/lib/python3.7/site-packages/pyscenic/utils.py", line 136, in add_correlation
    rhos = np.array([corr_mtx[s2][s1] for s1, s2 in zip(adjacencies.TF, adjacencies.target)])
  File "/home/Morgane/anaconda3/lib/python3.7/site-packages/pyscenic/utils.py", line 136, in <listcomp>
    rhos = np.array([corr_mtx[s2][s1] for s1, s2 in zip(adjacencies.TF, adjacencies.target)])
  File "/home/Morgane/anaconda3/lib/python3.7/site-packages/pandas/core/frame.py", line 2980, in __getitem__
    indexer = self.columns.get_loc(key)
  File "/home/Morgane/anaconda3/lib/python3.7/site-packages/pandas/core/indexes/base.py", line 2899, in get_loc
    return self._engine.get_loc(self._maybe_cast_indexer(key))
  File "pandas/_libs/index.pyx", line 107, in pandas._libs.index.IndexEngine.get_loc
  File "pandas/_libs/index.pyx", line 128, in pandas._libs.index.IndexEngine.get_loc
  File "pandas/_libs/index_class_helper.pxi", line 91, in pandas._libs.index.Int64Engine._check_type
KeyError: 'RPS19'

I have read other users errors and tried trouble shooting without any luck.

thanks for your help
Morgane 

memory error on SLURM servers
I have been using pySCENIC on a few datasets locally but have been looking into larger datasets that required me to move to an HPC infrastructure.
My dataset is not huge, around 10k cells and 17k genes. I have been able to successfully run the `grnboost2` step, and have saved the output for that.
However, I keep running into issues when I try pruning the modules. I have an `OSError: [Errno 12] Cannot allocate memory` error, as I show below, even though the job does not exceed 8% of the available memory usage.

1. Have you experienced this before, and how should I deal with this?
2. How can I derive regulons without pruning?

```
>>> df = prune2df(dbs, modules, MOTIF_ANNOTATIONS_FNAME)
Traceback (most recent call last):
  File "<stdin>", line 1, in <module>
  File "/accounts/campus/koenvdberge/.conda/envs/scenic2/lib/python3.7/site-packages/pyscenic/prune.py", line 351, in prune2df
    num_workers, module_chunksize)
  File "/accounts/campus/koenvdberge/.conda/envs/scenic2/lib/python3.7/site-packages/pyscenic/prune.py", line 300, in _distributed_calc
    return create_graph().compute(scheduler='processes', num_workers=num_workers if num_workers else cpu_count())
  File "/accounts/campus/koenvdberge/.conda/envs/scenic2/lib/python3.7/site-packages/dask/base.py", line 156, in compute
    (result,) = compute(self, traverse=False, **kwargs)
  File "/accounts/campus/koenvdberge/.conda/envs/scenic2/lib/python3.7/site-packages/dask/base.py", line 397, in compute
    results = schedule(dsk, keys, **kwargs)
  File "/accounts/campus/koenvdberge/.conda/envs/scenic2/lib/python3.7/site-packages/dask/multiprocessing.py", line 167, in get
    initializer=initialize_worker_process)
  File "/accounts/campus/koenvdberge/.conda/envs/scenic2/lib/python3.7/multiprocessing/context.py", line 119, in Pool
    context=self.get_context())
  File "/accounts/campus/koenvdberge/.conda/envs/scenic2/lib/python3.7/multiprocessing/pool.py", line 176, in __init__
    self._repopulate_pool()
  File "/accounts/campus/koenvdberge/.conda/envs/scenic2/lib/python3.7/multiprocessing/pool.py", line 241, in _repopulate_pool
    w.start()
  File "/accounts/campus/koenvdberge/.conda/envs/scenic2/lib/python3.7/multiprocessing/process.py", line 112, in start
    self._popen = self._Popen(self)
  File "/accounts/campus/koenvdberge/.conda/envs/scenic2/lib/python3.7/multiprocessing/context.py", line 277, in _Popen
    return Popen(process_obj)
  File "/accounts/campus/koenvdberge/.conda/envs/scenic2/lib/python3.7/multiprocessing/popen_fork.py", line 20, in __init__
    self._launch(process_obj)
  File "/accounts/campus/koenvdberge/.conda/envs/scenic2/lib/python3.7/multiprocessing/popen_fork.py", line 70, in _launch
    self.pid = os.fork()
OSError: [Errno 12] Cannot allocate memory
>>> dbs
[FeatherRankingDatabase(name="mm9-tss-centered-10kb-7species.mc9nr"), FeatherRankingDatabase(name="mm9-500bp-upstream-7species.mc9nr")]
>>> len(modules)
6801
>>> MOTIF_ANNOTATIONS_FNAME
'/accounts/campus/koenvdberge/Documents/SCENIC/motif2tf/mouse/motifs-v9-nr.mgi-m0.001-o0.0.tbl'
```


Discussion / hints:


Hi,

This issue looks similar to #103 , but I'm not sure if it was every really solved. Perhaps you could try running this step via the CLI?
Hello, 

Do you have a tutorial for the CLI? Can I start with the expression matrix filtered fro the R package or do I need to run correlation first and then export it?

thanks for your answer
Morgane

From [this tutorial](http://htmlpreview.github.io/?https://github.com/aertslab/SCENICprotocol/blob/master/notebooks/PBMC10k_SCENIC-protocol-CLI.html) (search for "cisTarget"), you can see the usage of the command line:

```
pyscenic ctx adj.tsv \
    {f_db_names} \
    --annotations_fname {f_motif_path} \
    --expression_mtx_fname {f_loom_path_scenic} \
    --output reg.csv \
    --num_workers 20
```

The correlation step is baked into this already, but you may need to write out your grn output (adj.tsv) to a text file if you're running from Jupyter.
Thanks for your help, 
I had to re run the first part to get the adjacencies file in the correct form. I started the second part to get regulons 2 days ago and it is still running. do you have an idea of how long it might take? my matrix is about 40,000 cells and more than 22,000 genes.
Morgane
@morganee261 , sorry for the delayed reply. It sounds like you fixed this particular problem so I'll close this
I have not used pySCENIC on slurm, but this must be something specific to your HPC, I would guess. It's clearly exceeding some memory limitation (a per process limit perhaps?). For this dataset, I'd expect maybe 2 or 3GB used per process, which isn't much. Have you tried lowering the number of processes? 

You could also try the command line version instead of running it interactively (from a shell, check out `pyscenic ctx -h`).

If you just want to look at the TF-gene modules from the GRN without running the pruning step, you should be able to get it from the previous command:
```
modules = list(modules_from_adjacencies(adjacencies, exprMat))
```
but this won't provide any of the benefit of module pruning to generate regulons, of course.
Thanks for the swift reply. I deliberately did not use any parallelization here in an attempt to try and find out what was going wrong, and since the memory tracking is showing me that memory consumption is actually rather low, I was thinking the error may originate from something else, rather than a true memory issue. I will look into the command line version.  
I am trying the CLI right now, where it starts running, but I seem to be getting an error related to auc. Any ideas?

```
(scenic2) gandalf.koenvdberge$ pyscenic ctx -o "/accounts/campus/koenvdberge/Documents/CRISPRi/pyScenic_cli/output2_prune" --annotations_fname "/accounts/campus/koenvdberge/Documents/SCENIC/motif2tf/mouse/motifs-v9-nr.mgi-m0.001-o0.0.tbl" "/accounts/campus/koenvdberge/Documents/CRISPRi/pyScenicResults/modules.dat" "/accounts/campus/koenvdberge/Documents/SCENIC/cisTarget_databases/mm9-tss-centered-10kb-7species.mc9nr.feather" "/accounts/campus/koenvdberge/Documents/SCENIC/cisTarget_databases/mm9-500bp-upstream-7species.mc9nr.feather"  

/accounts/campus/koenvdberge/.conda/envs/scenic2/lib/python3.7/site-packages/dask/config.py:161: YAMLLoadWarning: calling yaml.load() without Loader=... is deprecated, as the default Loader is unsafe. Please read https://msg.pyyaml.org/load for full details.
  data = yaml.load(f.read()) or {}

2020-04-20 20:48:59,451 - pyscenic.cli.pyscenic - INFO - Loading modules.

2020-04-20 20:49:01,867 - pyscenic.cli.pyscenic - INFO - Loading databases.

2020-04-20 20:49:01,867 - pyscenic.cli.pyscenic - INFO - Calculating regulons.
[                                        ] | 0% Completed | 25.1s
2020-04-20 20:49:29,803 - pyscenic.transform - WARNING - Less than 80% of the genes in Regulon for EEF1D could be mapped to mm9-tss-centered-10kb-7species.mc9nr. Skipping this module.
[                                        ] | 0% Completed | 27.3s
2020-04-20 20:49:31,903 - pyscenic.transform - WARNING - Less than 80% of the genes in Regulon for EGR1 could be mapped to mm9-tss-centered-10kb-7species.mc9nr. Skipping this module.
[                                        ] | 0% Completed | 39.8s
2020-04-20 20:49:44,471 - pyscenic.transform - WARNING - Less than 80% of the genes in Regulon for HNRNPC could be mapped to mm9-tss-centered-10kb-7species.mc9nr. Skipping this module.
[                                        ] | 0% Completed | 42.5s
2020-04-20 20:49:47,100 - pyscenic.transform - WARNING - Less than 80% of the genes in Regulon for HNRNPH3 could be mapped to mm9-tss-centered-10kb-7species.mc9nr. Skipping this module.
[                                        ] | 0% Completed | 47.4s
2020-04-20 20:49:52,063 - pyscenic.transform - WARNING - Less than 80% of the genes in Regulon for LARP1 could be mapped to mm9-tss-centered-10kb-7species.mc9nr. Skipping this module.
[                                        ] | 0% Completed | 51.8s
Traceback (most recent call last):
  File "/accounts/campus/koenvdberge/.conda/envs/scenic2/bin/pyscenic", line 8, in <module>
    sys.exit(main())
  File "/accounts/campus/koenvdberge/.conda/envs/scenic2/lib/python3.7/site-packages/pyscenic/cli/pyscenic.py", line 408, in main
    args.func(args)
  File "/accounts/campus/koenvdberge/.conda/envs/scenic2/lib/python3.7/site-packages/pyscenic/cli/pyscenic.py", line 159, in prune_targets_command
    num_workers=args.num_workers)
  File "/accounts/campus/koenvdberge/.conda/envs/scenic2/lib/python3.7/site-packages/pyscenic/prune.py", line 351, in prune2df
    num_workers, module_chunksize)
  File "/accounts/campus/koenvdberge/.conda/envs/scenic2/lib/python3.7/site-packages/pyscenic/prune.py", line 300, in _distributed_calc
    return create_graph().compute(scheduler='processes', num_workers=num_workers if num_workers else cpu_count())
  File "/accounts/campus/koenvdberge/.conda/envs/scenic2/lib/python3.7/site-packages/dask/base.py", line 156, in compute
    (result,) = compute(self, traverse=False, **kwargs)
  File "/accounts/campus/koenvdberge/.conda/envs/scenic2/lib/python3.7/site-packages/dask/base.py", line 397, in compute
    results = schedule(dsk, keys, **kwargs)
  File "/accounts/campus/koenvdberge/.conda/envs/scenic2/lib/python3.7/site-packages/dask/multiprocessing.py", line 192, in get
    raise_exception=reraise, **kwargs)
  File "/accounts/campus/koenvdberge/.conda/envs/scenic2/lib/python3.7/site-packages/dask/local.py", line 501, in get_async
    raise_exception(exc, tb)
  File "/accounts/campus/koenvdberge/.conda/envs/scenic2/lib/python3.7/site-packages/dask/compatibility.py", line 112, in reraise
    raise exc
  File "/accounts/campus/koenvdberge/.conda/envs/scenic2/lib/python3.7/site-packages/dask/local.py", line 272, in execute_task
    result = _execute_task(task, data)
  File "/accounts/campus/koenvdberge/.conda/envs/scenic2/lib/python3.7/site-packages/dask/local.py", line 252, in _execute_task
    args2 = [_execute_task(a, cache) for a in args]
  File "/accounts/campus/koenvdberge/.conda/envs/scenic2/lib/python3.7/site-packages/dask/local.py", line 252, in <listcomp>
    args2 = [_execute_task(a, cache) for a in args]
  File "/accounts/campus/koenvdberge/.conda/envs/scenic2/lib/python3.7/site-packages/dask/local.py", line 253, in _execute_task
    return func(*args2)
  File "/accounts/campus/koenvdberge/.conda/envs/scenic2/lib/python3.7/site-packages/pyscenic/transform.py", line 231, in modules2df
    for module in modules])
  File "/accounts/campus/koenvdberge/.conda/envs/scenic2/lib/python3.7/site-packages/pyscenic/transform.py", line 231, in <listcomp>
    for module in modules])
  File "/accounts/campus/koenvdberge/.conda/envs/scenic2/lib/python3.7/site-packages/pyscenic/transform.py", line 185, in module2df
    weighted_recovery=weighted_recovery)
  File "/accounts/campus/koenvdberge/.conda/envs/scenic2/lib/python3.7/site-packages/pyscenic/transform.py", line 129, in module2features_auc1st_impl
    aucs = calc_aucs(df, db.total_genes, weights, auc_threshold)
  File "/accounts/campus/koenvdberge/.conda/envs/scenic2/lib/python3.7/site-packages/pyscenic/recovery.py", line 284, in aucs
    assert maxauc > 0
AssertionError
```
It looks like your gene matrix contains human (HGNC gene symbols) gene names instead of mouse (MGI gene symbols), so in that case use the human databases.

Hence the `Less than 80% of the genes in Regulon for EEF1D could be mapped to mm9-tss-centered-10kb-7species.mc9nr` warnings.
Thanks, working with the command line interface instead worked for me in avoiding the errors.

You were right that I had selected the mouse database while I should have selected the human one. However, a lot of these  `Less than 80% of the genes` messages keep popping up which results in a low number of TFs after pruning. Is there a way TFs must be defined, eg uppercase vs lowercase or something I could check?
Gene names for Human should be in HGNC gene symbols (uppercase):
- https://www.genenames.org
- https://biomart.genenames.org/martform/#!/default/HGNC?datasets=hgnc_gene_mart&upload_list%5B%5D=hgnc_gene__hgnc_gene_id_1010&upload_list%5B%5D=&upload_list%5B%5D=HGNC+ID(s)
Thanks, that's indeed what I am using, but I am still getting abundant warnings. Is that expected?

Here is my code

```
pyscenic grn -t -m grnboost2 -o "/accounts/campus/koenvdberge/Documents/CRISPRi/pyScenic_cli/output1_grn" "/accounts/campus/koenvdberge/Documents/CRISPRi/data/counts_nogRNA.csv" "/accounts/campus/koenvdberge/Documents/CRISPRi/tf_names.csv"

pyscenic ctx -o "/accounts/campus/koenvdberge/Documents/CRISPRi/pyScenic_cli/output2_prune.dat" --annotations_fname "/accounts/campus/koenvdberge/Documents/SCENIC/motif2tf/human/motifs-v9-nr.hgnc-m0.001-o0.0.tbl" "/accounts/campus/koenvdberge/Documents/CRISPRi/pyScenic_cli/output1_grn.csv" "/accounts/campus/koenvdberge/Documents/SCENIC/cisTarget_databases/human/hg19-500bp-upstream-7species.mc9nr.feather" "/accounts/campus/koenvdberge/Documents/SCENIC/cisTarget_databases/human/hg19-tss-centered-10kb-7species.mc9nr.feather"

pyscenic aucell -t -o "/accounts/campus/koenvdberge/Documents/CRISPRi/pyScenic_cli/output3_aucell" "/accounts/campus/koenvdberge/Documents/CRISPRi/data/counts_nogRNA.csv" "/accounts/campus/koenvdberge/Documents/CRISPRi/pyScenic_cli/output2_prune.dat"
```
Could you post a section of your gene names?
Sure, here's some gene names:

```
 "TTTY17B"   "TTTY17C"   "TTTY17A"   "TTTY4C"    "TTTY4B"    "TTTY4"     "BPY2C"     "BPY2B"     "BPY2" "DAZ4"      "DAZ1"      "DAZ3"      "DAZ2"      "TTTY3B"    "TTTY3"     "CDY1B"     "CDY1"      "CSPG4P1Y" "GOLGA2P3Y" "GOLGA2P2Y"
```

and TF names:

```
HOXA9
ZFP128
ZFP853
NR1H2
NR1H3
NR1H4
NR1H5
NR1I2
```
I think the issue stems from the fact that you have long non-coding RNAs in your expression matrix. Our databases only contain normal genes (mostly protein coding, there might be some pseudogenes in there too), but not lincRNAs, microRNAs, ...RNA.
https://www.genenames.org/tools/search/#!/all?query=TTTY4C
Thanks, so I suspect it would be better to remove these prior to running `pySCENIC` such that the regulons aren't 'contaminated' with non-protein coding genes.
Yes remove them before running the first step.

Fix patch:
diff --git a/.travis.yml b/.travis.yml
index 79504b0..942b734 100644
--- a/.travis.yml
+++ b/.travis.yml
@@ -1,12 +1,14 @@
+dist: bionic
+
 language: python
-python:
-    - "3.6"
+
+matrix:
+  include:
+    - python: 3.6
+    - python: 3.7
+    - python: 3.8
+
 install: pip install tox-travis
+
 script: tox
-#deploy:
-#  provider: pypi
-#  user: ""
-#  password:
-#    secure: ""
-#  on:
-#    tags: true
\ No newline at end of file
+
diff --git a/README.rst b/README.rst
index e206e06..483a08a 100644
--- a/README.rst
+++ b/README.rst
@@ -18,30 +18,49 @@ in no time. The latter is achieved via the dask_ framework for distributed compu
 News and releases
 -----------------
 
-0.10.4 | 2020-11-24
-^^^^^^^^^^^^^^^^^^^
+0.11.0 | 2021-02-10
+^^^^^^^^^^^^^^^^
 
-* Included new CLI option to add correlation information to the GRN adjacencies file. This can be called with ``pyscenic add_cor``.
+**Major features:**
 
-0.10.3 | 2020-07-15
-^^^^^^^^^^^^^^^^^^^
+* Updated `Arboreto <https://github.com/aertslab/arboreto>`_ release (GRN inference step) includes:
 
-* Integrate arboreto multiprocessing script into pySCENIC CLI
-* Skip modules with zero db overlap in cisTarget step
-* Additional error message if regulons file is empty
-* Additional error if there is a mismatch between the genes present in the GRN and the expression matrix
-* Fixed bug in motif url construciton when running without pruning
+  * Support for sparse matrices (using the ``--sparse`` flag in ``pyscenic grn``, or passing a sparse matrix to ``grnboost2``/``genie3``).
+  * Fixes to avoid dask metadata mismatch error
 
-0.10.2 | 2020-06-05
-^^^^^^^^^^^^^^^^^^^
+* Updated cisTarget:
+
+  * Fix for metadata mismatch in ctx prune2df step
+  * Support for databases Apache Parquet format
+  * Faster loading from feather databases
+  * Bugfix: loading genes from a database (previously missing the last gene name in the database)
+
+* Support for Anndata input and output
+
+* Package updates:
 
-* Bugfix for CLI grn step
+  * Upgrade to newer pandas version
+  * Upgrade to newer numba version
+  * Upgrade to newer versions of dask, distributed
 
+* Input checks and more descriptive error messages.
 
-0.10.1 | 2020-05-17
+  * Check that regulons loaded are not empty.
+
+* Bugfixes:
+
+  * Motif url construction fixed when running ctx without pruning
+  * Compression of intermediate files in the CLI steps
+  * Handle loom files with non-standard gene/cell attribute names
+  * Reformat the genesig gmt input/output
+  * Fix AUCell output to loom with non-standard loom attributes
+
+
+0.10.4 | 2020-11-24
 ^^^^^^^^^^^^^^^^^^^
 
-* CLI: file compression (optionally) enabled for intermediate files for the major steps: grn (adjacencies matrix), ctx (regulons), and aucell (auc matrix). Compression is used when the file name argument has a .gz ending.
+* Included new CLI option to add correlation information to the GRN adjacencies file. This can be called with ``pyscenic add_cor``.
+
 
 
 See also the extended `Release Notes <https://pyscenic.readthedocs.io/en/latest/releasenotes.html>`_.
@@ -94,8 +113,8 @@ References
 .. |buildstatus| image:: https://travis-ci.org/aertslab/pySCENIC.svg?branch=master
 .. _buildstatus: https://travis-ci.org/aertslab/pySCENIC
 
-.. |pypipackage| image:: https://badge.fury.io/py/pyscenic.svg
-.. _pypipackage: https://badge.fury.io/py/pyscenic
+.. |pypipackage| image:: https://img.shields.io/pypi/v/pySCENIC?color=%23026aab
+.. _pypipackage: https://pypi.org/project/pyscenic/
 
 .. |docstatus| image:: https://readthedocs.org/projects/pyscenic/badge/?version=latest
 .. _docstatus: http://pyscenic.readthedocs.io/en/latest/?badge=latest
diff --git a/docs/releasenotes.rst b/docs/releasenotes.rst
index d4c6c45..8f40837 100644
--- a/docs/releasenotes.rst
+++ b/docs/releasenotes.rst
@@ -1,6 +1,29 @@
 Release Notes
 =============
 
+
+0.10.3 | 2020-07-15
+^^^^^^^^^^^^^^^^^^^
+
+* Integrate arboreto multiprocessing script into pySCENIC CLI
+* Skip modules with zero db overlap in cisTarget step
+* Additional error message if regulons file is empty
+* Additional error if there is a mismatch between the genes present in the GRN and the expression matrix
+* Fixed bug in motif url construction when running without pruning
+
+
+0.10.2 | 2020-06-05
+^^^^^^^^^^^^^^^^^^^
+
+* Bugfix for CLI grn step
+
+
+0.10.1 | 2020-05-17
+^^^^^^^^^^^^^^^^^^^
+
+* CLI: file compression (optionally) enabled for intermediate files for the major steps: grn (adjacencies matrix), ctx (regulons), and aucell (auc matrix). Compression is used when the file name argument has a .gz ending.
+
+
 0.10.3 | 2020-07-15
 ^^^^^^^^^^^^^^^^^^^
 
@@ -10,6 +33,7 @@ Release Notes
 * Additional error if there is a mismatch between the genes present in the GRN and the expression matrix
 * Fixed bug in motif url construciton when running without pruning
 
+
 0.10.2 | 2020-06-05
 ^^^^^^^^^^^^^^^^^^^
 
diff --git a/requirements.txt b/requirements.txt
index d931e4d..73c49af 100644
--- a/requirements.txt
+++ b/requirements.txt
@@ -1,16 +1,16 @@
 cytoolz
 multiprocessing_on_dill
 llvmlite
-numba
+numba>=0.51.2
 attrs
 frozendict
 numpy
-pandas>=0.20.1,<1.0.0
+pandas>=0.20.1
 cloudpickle
-dask==1.0.0
-distributed>=1.21.6,<2.0.0
+dask
+distributed
 pyarrow>=0.11.1,<0.17.0
-arboreto
+arboreto>=0.1.6
 boltons
 setuptools
 pyyaml
@@ -20,3 +20,6 @@ umap-learn
 loompy
 networkx
 scipy
+fsspec
+requests
+aiohttp
diff --git a/requirements_docker.txt b/requirements_docker.txt
index 801a9c6..43ff6f5 100644
--- a/requirements_docker.txt
+++ b/requirements_docker.txt
@@ -60,7 +60,7 @@ networkx==2.3
 numba==0.46.0
 numexpr==2.7.0
 numpy==1.17.2
-pandas==0.25.3
+pandas==1.0.3
 pandocfilters==1.4.2
 papermill==1.2.0
 parso==0.5.1
diff --git a/src/pyscenic/cli/db2feather.py b/src/pyscenic/cli/db2feather.py
index 68f8cd7..9e80109 100644
--- a/src/pyscenic/cli/db2feather.py
+++ b/src/pyscenic/cli/db2feather.py
@@ -2,7 +2,7 @@
 
 import os
 import argparse
-from pyscenic.rnkdb import convert2feather
+from pyscenic.rnkdb import convert_sqlitedb_to_featherdb
 
 
 def derive_db_name(fname:str) -> str:
@@ -25,7 +25,7 @@ def create_argument_parser():
 def convert(out_folder, in_fnames):
     for fname in in_fnames:
         print("Converting {}".format(fname.name))
-        convert2feather(fname.name, out_folder, derive_db_name(fname.name))
+        convert_sqlitedb_to_featherdb(fname.name, out_folder, derive_db_name(fname.name))
 
 
 def main():
diff --git a/src/pyscenic/cli/pyscenic.py b/src/pyscenic/cli/pyscenic.py
index 5cfcbf4..1da8e87 100644
--- a/src/pyscenic/cli/pyscenic.py
+++ b/src/pyscenic/cli/pyscenic.py
@@ -71,7 +71,10 @@ def find_adjacencies_command(args):
     client, shutdown_callback = _prepare_client(args.client_or_address, num_workers=args.num_workers)
     method = grnboost2 if args.method == 'grnboost2' else genie3
     try:
-        network = method(expression_data=ex_mtx, tf_names=tf_names, verbose=True, client_or_address=client, seed=args.seed)
+        if args.sparse:
+            network = method(expression_data=ex_mtx[0], gene_names=ex_mtx[1], tf_names=tf_names, verbose=True, client_or_address=client, seed=args.seed)
+        else:
+            network = method(expression_data=ex_mtx, tf_names=tf_names, verbose=True, client_or_address=client, seed=args.seed)
     finally:
         shutdown_callback(False)
 
@@ -234,10 +237,22 @@ def aucell_command(args):
     if '.loom' in extension:
         try:
             copyfile(args.expression_mtx_fname.name, args.output.name)
-            append_auc_mtx(args.output.name, auc_mtx, signatures, args.seed, args.num_workers)
+            append_auc_mtx(args.output.name, ex_mtx, auc_mtx, signatures, args.seed, args.num_workers)
         except OSError as e:
             LOGGER.error("Expression matrix should be provided in the loom file format.")
             sys.exit(1)
+    elif '.h5ad' in extension:
+        from pyscenic.export import add_scenic_metadata
+        from anndata import read_h5ad
+        # check input file is also h5ad:
+        if '.h5ad' in PurePath(args.expression_mtx_fname.name).suffixes:
+            copyfile(args.expression_mtx_fname.name, args.output.name)
+            add_scenic_metadata(read_h5ad(filename=args.output.name, backed='r'),
+                                auc_mtx,
+                                signatures).write(args.output.name)
+        else:
+            LOGGER.error("Expression matrix should be provided in the h5ad (anndata) file format.")
+            sys.exit(1)
     elif args.output.name == '<stdout>':
         transpose = (args.transpose == 'yes')
         (auc_mtx.T if transpose else auc_mtx).to_csv(args.output)
diff --git a/src/pyscenic/cli/utils.py b/src/pyscenic/cli/utils.py
index 4928776..5908158 100644
--- a/src/pyscenic/cli/utils.py
+++ b/src/pyscenic/cli/utils.py
@@ -85,7 +85,7 @@ def suffixes_to_separator(extension):
 def is_valid_suffix(extension, method):
     assert(isinstance(extension,list)), 'extension should be of type "list"'
     if method in ['grn', 'aucell']:
-        valid_extensions = ['.csv', '.tsv', '.loom']
+        valid_extensions = ['.csv', '.tsv', '.loom', '.h5ad']
     elif method == 'ctx':
         valid_extensions = ['.csv', '.tsv']
     elif method == 'ctx_yaml':
@@ -114,6 +114,15 @@ def load_exp_matrix(fname: str, transpose: bool = False,
     if is_valid_suffix(extension, 'grn'):
         if '.loom' in extension:
             return load_exp_matrix_as_loom(fname, return_sparse, attribute_name_cell_id, attribute_name_gene)
+        elif '.h5ad' in extension:
+            from anndata import read_h5ad
+            adata = read_h5ad(filename=fname, backed='r')
+            if return_sparse:
+                # expr, gene, cell:
+                return adata.X.value, adata.var_names.values, adata.obs_names.values
+            else:
+                return pd.DataFrame(adata.X.value.todense(), index=adata.obs_names.values, columns=adata.var_names.values)
+
         else:
             df = pd.read_csv(fname, sep=suffixes_to_separator(extension), header=0, index_col=0)
             return df.T if transpose else df
@@ -180,7 +189,7 @@ def load_signatures(fname: str) -> Sequence[Type[GeneSignature]]:
         return GeneSignature.from_gmt(fname,
                                   field_separator=sep,
                                   gene_separator=sep)
-    elif extension == '.dat':
+    elif '.dat' in extension:
         with openfile(fname, 'rb') as f:
             return pickle.load(f)
     else:
@@ -234,10 +243,7 @@ def load_modules(fname: str) -> Sequence[Type[GeneSignature]]:
         with openfile(fname, 'rb') as f:
             return pickle.load(f)
     elif '.gmt' in extension:
-        sep = guess_separator(fname)
-        return GeneSignature.from_gmt(fname,
-                                      field_separator=sep,
-                                      gene_separator=sep)
+        return GeneSignature.from_gmt(fname)
     else:
         raise ValueError("Unknown file format for \"{}\".".format(fname))
 
@@ -254,7 +260,7 @@ def compress_meta(meta):
     return base64.b64encode(zlib.compress(json.dumps(meta).encode('ascii'))).decode('ascii')
 
 
-def append_auc_mtx(fname: str, auc_mtx: pd.DataFrame, regulons: Sequence[Type[GeneSignature]], seed=None, num_workers=1) -> None:
+def append_auc_mtx(fname: str, ex_mtx: pd.DataFrame, auc_mtx: pd.DataFrame, regulons: Sequence[Type[GeneSignature]], seed=None, num_workers=1) -> None:
     """
 
     Append AUC matrix to loom file.
@@ -283,7 +289,6 @@ def fetch_logo(context):
                            "motifData": name2logo.get(name, "")} for name, threshold in auc_thresholds.iteritems()]
 
     # Calculate the number of genes per cell.
-    ex_mtx = load_exp_matrix(fname)
     binary_mtx = ex_mtx.copy()
     binary_mtx[binary_mtx != 0] = 1.0
     ngenes = binary_mtx.sum(axis=1).astype(int)
diff --git a/src/pyscenic/genesig.py b/src/pyscenic/genesig.py
index 815660d..305a212 100644
--- a/src/pyscenic/genesig.py
+++ b/src/pyscenic/genesig.py
@@ -61,7 +61,7 @@ def from_yaml(cls, loader, node):
                              gene2weight=zip(data['genes'], data['weights']))
 
     @classmethod
-    def from_gmt(cls, fname: str, field_separator: str = ',', gene_separator: str = ',') -> List['GeneSignature']:
+    def from_gmt(cls, fname: str, field_separator: str = '\t', gene_separator: str = '\t') -> List['GeneSignature']:
         """
         Load gene signatures from a GMT file.
 
@@ -81,13 +81,13 @@ def signatures():
                     if line.startswith("#") or not line.strip():
                         continue
                     columns = re.split(field_separator, line.rstrip())
-                    genes = columns[2:] if field_separator == gene_separator else columns[2].split(gene_separator)
+                    genes = columns[2:]
                     yield GeneSignature(name=columns[0], gene2weight=genes)
         return list(signatures())
 
 
     @classmethod
-    def to_gmt(cls, fname: str, signatures: List[Type['GeneSignature']], field_separator: str = ',', gene_separator: str = ',') -> None:
+    def to_gmt(cls, fname: str, signatures: List[Type['GeneSignature']], field_separator: str = '\t', gene_separator: str = '\t') -> None:
         """
         Save list of signatures as GMT file.
 
@@ -101,7 +101,7 @@ def to_gmt(cls, fname: str, signatures: List[Type['GeneSignature']], field_separ
             for signature in signatures:
                 genes = gene_separator.join(signature.genes)
                 file.write("{}{}{}{}{}\n".format(signature.name, field_separator,
-                                                 signature.metadata(gene_separator), field_separator,
+                                                 signature.metadata(','), field_separator,
                                                  genes))
 
 
diff --git a/src/pyscenic/math.py b/src/pyscenic/math.py
index 513ebba..5cceb02 100644
--- a/src/pyscenic/math.py
+++ b/src/pyscenic/math.py
@@ -2,10 +2,10 @@
 
 
 import numpy as np
-from numba import *
+from numba import njit, float64, int64, prange
 
 
-@njit(signature_or_function=float64(float64[:], float64[:], float64))
+@njit(float64(float64[:], float64[:], float64))
 def masked_rho(x: np.ndarray, y: np.ndarray, mask: float = 0.0) -> float:
     """
     Calculates the masked correlation coefficient of two vectors.
@@ -30,7 +30,7 @@ def masked_rho(x: np.ndarray, y: np.ndarray, mask: float = 0.0) -> float:
     return cov_xy / (std_x * std_y)
 
 
-@njit(signature_or_function=float64[:, :](float64[:, :], float64[:, :], float64), parallel=True)
+@njit(float64[:, :](float64[:, :], float64[:, :], float64), parallel=True)
 def masked_rho_2d(x: np.ndarray, y: np.ndarray, mask: float = 0.0) -> np.ndarray:
     """
     Calculates the masked correlation coefficients of two arrays.
@@ -50,7 +50,7 @@ def masked_rho_2d(x: np.ndarray, y: np.ndarray, mask: float = 0.0) -> np.ndarray
     return rhos
 
 
-@njit(signature_or_function=float64[:](float64[:, :], int64[:, :], float64), parallel=True)
+@njit(float64[:](float64[:, :], int64[:, :], float64), parallel=True)
 def masked_rho4pairs(mtx: np.ndarray, col_idx_pairs: np.ndarray, mask: float = 0.0) -> np.ndarray:
     """
     Calculates the masked correlation of columns pairs in a matrix.
@@ -66,4 +66,5 @@ def masked_rho4pairs(mtx: np.ndarray, col_idx_pairs: np.ndarray, mask: float = 0
         x = mtx[:, col_idx_pairs[n_idx, 0]]
         y = mtx[:, col_idx_pairs[n_idx, 1]]
         rhos[n_idx] = masked_rho(x, y, mask)
-    return rhos
\ No newline at end of file
+    return rhos
+
diff --git a/src/pyscenic/recovery.py b/src/pyscenic/recovery.py
index 8f2ab27..66eda38 100644
--- a/src/pyscenic/recovery.py
+++ b/src/pyscenic/recovery.py
@@ -4,7 +4,7 @@
 import numpy as np
 from itertools import repeat
 from typing import Type, Optional, List, Tuple
-from numba import *
+from numba import jit
 import logging
 
 from .rnkdb import RankingDatabase
diff --git a/src/pyscenic/rnkdb.py b/src/pyscenic/rnkdb.py
index 74fff69..0a616b9 100644
--- a/src/pyscenic/rnkdb.py
+++ b/src/pyscenic/rnkdb.py
@@ -1,17 +1,52 @@
 # -*- coding: utf-8 -*-
 
 import os
-import pandas as pd
 import numpy as np
-from typing import Tuple, Set, Type
-from abc import ABCMeta, abstractmethod
+import pandas as pd
+import pyarrow as pa
+import pyarrow.parquet as pq
 import sqlite3
+
+from abc import ABCMeta, abstractmethod
 from operator import itemgetter
-from .genesig import GeneSignature
+from typing import Dict, Set, Tuple, Type
+
 from cytoolz import memoize
 from pyarrow.feather import write_feather, FeatherReader
 from tqdm import tqdm
 
+from .genesig import GeneSignature
+
+
+class PyArrowThreads:
+    """
+    A static class to control how many threads PyArrow is allowed to use to convert a Feather database to a pandas
+    dataframe.
+
+    By default the number of threads is set to 4.
+    Overriding the number of threads is possible by using the environment variable "PYARROW_THREADS=nbr_threads".
+    """
+    pyarrow_threads = 4
+
+    if os.environ.get("PYARROW_THREADS"):
+        try:
+            # If "PYARROW_THREADS" is set, check if it is a number.
+            pyarrow_threads = int(os.environ.get("PYARROW_THREADS"))
+        except ValueError:
+            pass
+
+        if pyarrow_threads < 1:
+            # Set the number of PyArrow threads to 1 if a negative number or zero was specified.
+            pyarrow_threads = 1
+
+    @staticmethod
+    def set_nbr_pyarrow_threads(nbr_threads=None):
+        # Set number of threads to use for PyArrow when converting Feather database to pandas dataframe.
+        pa.set_cpu_count(nbr_threads if nbr_threads else PyArrowThreads.pyarrow_threads)
+
+
+PyArrowThreads.set_nbr_pyarrow_threads()
+
 
 class RankingDatabase(metaclass=ABCMeta):
     """
@@ -243,6 +278,7 @@ def __init__(self, fname: str, name: str):
     @property
     @memoize
     def total_genes(self) -> int:
+        # Do not count column 1 as it contains the index with the name of the features.
         return FeatherReader(self._fname).num_columns - 1
 
     @property
@@ -250,13 +286,96 @@ def total_genes(self) -> int:
     def genes(self) -> Tuple[str]:
         # noinspection PyTypeChecker
         reader = FeatherReader(self._fname)
-        return tuple(reader.get_column_name(idx) for idx in range(self.total_genes) if reader.get_column_name(idx) != INDEX_NAME)
+        # Get all gene names (exclude "features" column).
+        return tuple(
+            reader.get_column_name(idx) for idx in range(reader.num_columns)
+            if reader.get_column_name(idx) != INDEX_NAME
+        )
+
+    @property
+    @memoize
+    def genes2idx(self) -> Dict[str, int]:
+        return {gene: idx for idx, gene in enumerate(self.genes)}
+
+    def load_full(self) -> pd.DataFrame:
+        df = FeatherReader(self._fname).read_pandas()
+        # Avoid copying the whole dataframe by replacing the index in place.
+        # This makes loading a database twice as fast in case the database file is already in the filesystem cache.
+        df.set_index(INDEX_NAME, inplace=True)
+        return df
+
+    def load(self, gs: Type[GeneSignature]) -> pd.DataFrame:
+        # For some genes in the signature there might not be a rank available in the database.
+        gene_set = self.geneset.intersection(set(gs.genes))
+        # Read ranking columns for genes in order they appear in the Feather file.
+        df = FeatherReader(self._fname).read_pandas(
+            columns=(INDEX_NAME,) + tuple(sorted(gene_set, key=lambda gene: self.genes2idx[gene]))
+        )
+        # Avoid copying the whole dataframe by replacing the index in place.
+        # This makes loading a database twice as fast in case the database file is already in the filesystem cache.
+        df.set_index(INDEX_NAME, inplace=True)
+        return df
+
+
+class ParquetRankingDatabase(RankingDatabase):
+    def __init__(self, fname: str, name: str):
+        """
+        Create a new parquet database.
+
+        :param fname: The filename of the database.
+        :param name: The name of the database.
+        """
+        super().__init__(name=name)
+
+        assert os.path.isfile(fname), "Database {0:s} doesn't exist.".format(fname)
+        # FeatherReader cannot be pickle (important for dask framework) so filename is field instead.
+        self._fname = fname
+
+    @property
+    @memoize
+    def total_genes(self) -> int:
+        # Do not count column 1 as it contains the index with the name of the features.
+        return pq.read_metadata(self._fname).num_columns - 1
+
+    @property
+    @memoize
+    def genes(self) -> Tuple[str]:
+        # noinspection PyTypeChecker
+        metadata = pq.read_metadata(self._fname)
+        assert metadata.num_row_groups == 1, \
+            "Parquet database {0:s} has more than one row group.".format(self._fname)
+        metadata_row_group = metadata.row_group(0)
+        # Get all gene names (exclude "features" column).
+        return tuple(
+            metadata_row_group.column(idx).path_in_schema
+            for idx in range(0, metadata.num_columns)
+            if metadata_row_group.column(idx).path_in_schema != INDEX_NAME
+        )
+
+    @property
+    @memoize
+    def genes2idx(self) -> Dict[str, int]:
+        return {gene: idx for idx, gene in enumerate(self.genes)}
 
     def load_full(self) -> pd.DataFrame:
-        return FeatherReader(self._fname).read_pandas().set_index(INDEX_NAME)
+        df = pq.read_pandas(self._fname).to_pandas()
+        # Avoid copying the whole dataframe by replacing the index in place.
+        # This makes loading a database twice as fast in case the database file is already in the filesystem cache.
+        df.set_index(INDEX_NAME, inplace=True)
+        return df
 
     def load(self, gs: Type[GeneSignature]) -> pd.DataFrame:
-        return FeatherReader(self._fname).read_pandas(columns=(INDEX_NAME,) + gs.genes).set_index(INDEX_NAME)
+        # For some genes in the signature there might not be a rank available in the database.
+        gene_set = self.geneset.intersection(set(gs.genes))
+        # Read rankings columns for genes in order they appear in the Parquet file.
+        df = pq.read_pandas(
+            self._fname,
+            columns=(INDEX_NAME,) + tuple(sorted(gene_set, key=lambda gene: self.genes2idx[gene]))
+        ).to_pandas()
+        # Avoid copying the whole dataframe by replacing the index in place.
+        # This makes loading a database twice as fast in case the database file is already in the filesystem cache.
+        df.set_index(INDEX_NAME, inplace=True)
+        return df
 
 
 class MemoryDecorator(RankingDatabase):
@@ -396,7 +515,10 @@ def __init__(self, fname: str, name: str):
         self.idx2identifier = {idx: identifier for identifier, idx in self.identifier2idx.items()}
 
         # Load dataframe into memory in a format most suited for fast loading of gene signatures.
-        df = FeatherReader(fname).read().set_index(INDEX_NAME)
+        df = FeatherReader(fname).read_pandas()
+        # Avoid copying the whole dataframe by replacing the index in place.
+        # This makes loading a database twice as fast in case the database file is already in the filesystem cache.
+        df.set_index(INDEX_NAME, inplace=True)
         self.max_rank = len(df.columns)
         self.features = [pd.Series(index=row.values, data=row.index, name=name) for name, row in df.iterrows()]
 
@@ -430,14 +552,14 @@ def load(self, gs: Type[GeneSignature]) -> pd.DataFrame:
                          axis=1).T.astype(INVERTED_DB_DTYPE).rename(columns=self.idx2identifier)
 
 
-def convert2feather(fname: str, out_folder: str, name: str, extension: str="feather") -> str:
+def convert_sqlitedb_to_featherdb(fname: str, out_folder: str, name: str, extension: str= "feather") -> str:
     """
-    Convert a whole genome rankings database to a feather format based database.
+    Convert a whole genome SQLite rankings database to a feather format based database.
 
     More information on this format can be found here:
     .. feather-format: https://blog.rstudio.com/2016/03/29/feather/
 
-    :param fname: The filename of the legacy
+    :param fname: The filename of the legacy SQLite rankings database.
     :param out_folder: The name of the folder to write the new database to.
     :param name: The name of the rankings database.
     :param extension: The extension of the new database file.
@@ -454,7 +576,8 @@ def convert2feather(fname: str, out_folder: str, name: str, extension: str="feat
     db = SQLiteRankingDatabase(fname=fname, name=name)
     df = db.load_full()
     df.index.name = INDEX_NAME
-    df.reset_index(inplace=True) # Index is not stored in feather format. https://github.com/wesm/feather/issues/200
+    # Index is not stored in feather format: https://github.com/wesm/feather/issues/200
+    df.reset_index(inplace=True)
     write_feather(df, feather_fname)
     return feather_fname
 
@@ -471,7 +594,10 @@ def opendb(fname: str, name: str) -> Type['RankingDatabase']:
     assert name, "A database should be given a proper name."
 
     extension = os.path.splitext(fname)[1]
-    if extension == ".feather":
+    if extension == ".parquet":
+        # noinspection PyTypeChecker
+        return ParquetRankingDatabase(fname, name=name)
+    elif extension == ".feather":
         if fname.endswith(".inverted.feather"):
             # noinspection PyTypeChecker
             return InvertedRankingDatabase(fname, name=name)
diff --git a/src/pyscenic/transform.py b/src/pyscenic/transform.py
index acdf964..851ef5c 100644
--- a/src/pyscenic/transform.py
+++ b/src/pyscenic/transform.py
@@ -224,7 +224,10 @@ def module2df(db: Type[RankingDatabase], module: Regulon, motif_annotations: pd.
     del df['Ranking']
     if not return_recovery_curves:
         del df['Recovery']
-    return df
+        assert all([ col in df.columns for col in DF_META_DATA ]), f"Column comparison to expected metadata failed! Found:\n{df.columns}"
+        return df[DF_META_DATA.columns]
+    else:
+        return df
 
 
 def modules2df(db: Type[RankingDatabase], modules: Sequence[Regulon], motif_annotations: pd.DataFrame,
diff --git a/src/pyscenic/utils.py b/src/pyscenic/utils.py
index 6fb9784..74528c7 100644
--- a/src/pyscenic/utils.py
+++ b/src/pyscenic/utils.py
@@ -323,7 +323,7 @@ def add_motif_url(df: pd.DataFrame, base_url: str):
     :param base_url:
     :return:
     """
-    df[("Enrichment", COLUMN_NAME_MOTIF_URL)] = list(map(partial(urljoin, base_url), df.index.get_level_values(COLUMN_NAME_MOTIF_ID)))
+    df[("Enrichment", COLUMN_NAME_MOTIF_URL)] = list(map(partial(urljoin, base=base_url), df.index.get_level_values(COLUMN_NAME_MOTIF_ID)))
     return df
 
 
diff --git a/tox.ini b/tox.ini
index a837a56..5ff038b 100644
--- a/tox.ini
+++ b/tox.ini
@@ -4,7 +4,7 @@
 # and then run "tox" from this directory.
 
 [tox]
-envlist = py36
+envlist = py36,py37,py38
 
 [testenv]
 commands = pytest


Test patch:
diff --git a/tests/test_aucell.py b/tests/test_aucell.py
index e37531d..aa2a230 100644
--- a/tests/test_aucell.py
+++ b/tests/test_aucell.py
@@ -8,8 +8,6 @@
 from pyscenic.aucell import derive_auc_threshold, aucell, create_rankings
 from pkg_resources import resource_filename
 
-
-
 NOMENCLATURE = "HGNC"
 TEST_EXPRESSION_MTX_FNAME = resource_filename('resources.tests', "GSE103322.em.hgnc.sample.cxg.csv")
 TEST_SIGNATURE_FNAME = resource_filename('resources.tests', "c6.all.v6.1.symbols.gmt")
@@ -22,15 +20,14 @@ def exp_matrix():
 
 @pytest.fixture
 def gs():
-    return GeneSignature.from_gmt(TEST_SIGNATURE_FNAME,
-                                  gene_separator="\t", field_separator="\t", )
+    return GeneSignature.from_gmt(TEST_SIGNATURE_FNAME, gene_separator="\t", field_separator="\t", )
 
 
 def test_create_rankings(exp_matrix):
     df_rnk = create_rankings(exp_matrix)
     n_genes = exp_matrix.shape[1]
     assert len(df_rnk.sum(axis=1).unique()) == 1
-    assert (df_rnk + 1).sum(axis=1).unique()[0] == (n_genes * (n_genes+1))/2.0
+    assert (df_rnk + 1).sum(axis=1).unique()[0] == (n_genes * (n_genes + 1)) / 2.0
 
 
 def test_aucell_w1(exp_matrix, gs):
@@ -48,6 +45,3 @@ def test_aucell_mismatch(exp_matrix, gs):
     gss = [GeneSignature(name="test", gene2weight=list(map("FAKE{}".format, range(100))))] + gs
     aucs_mtx = aucell(exp_matrix, gss, auc_threshold=percentiles[0.01], num_workers=1)
     print(aucs_mtx.head())
-
-
-
diff --git a/tests/test_featherdb.py b/tests/test_featherdb.py
index 8c5a64b..40eebf0 100644
--- a/tests/test_featherdb.py
+++ b/tests/test_featherdb.py
@@ -1,36 +1,44 @@
 # -*- coding: utf-8 -*-
 
 import pytest
-from pyscenic.rnkdb import FeatherRankingDatabase as RankingDatabase
-from pyscenic.genesig import GeneSignature
 from pkg_resources import resource_filename
 
+from pyscenic.genesig import GeneSignature
+from pyscenic.rnkdb import FeatherRankingDatabase as RankingDatabase
 
 TEST_DATABASE_FNAME = resource_filename('resources.tests', "hg19-tss-centered-10kb-10species.mc9nr.feather")
 TEST_DATABASE_NAME = "hg19-tss-centered-10kb-10species.mc9nr"
-TEST_SIGNATURE_FNAME = resource_filename('resources.tests', "/c6.all.v6.1.symbols.gmt")
+TEST_SIGNATURE_FNAME = resource_filename('resources.tests', "c6.all.v6.1.symbols.gmt")
 
 
 @pytest.fixture
 def db():
     return RankingDatabase(TEST_DATABASE_FNAME, TEST_DATABASE_NAME)
 
+
 @pytest.fixture
 def gs():
-    return GeneSignature.from_gmt(TEST_SIGNATURE_FNAME,
-                                    gene_separator="\t", field_separator="\t", )[0]
+    return GeneSignature.from_gmt(TEST_SIGNATURE_FNAME, gene_separator="\t", field_separator="\t", )[0]
+
 
 def test_init(db):
     assert db.name == TEST_DATABASE_NAME
 
+
 def test_total_genes(db):
     assert db.total_genes == 22284
 
+
+def test_genes(db):
+    assert len(db.genes) == 22284
+
+
 def test_load_full(db):
     rankings = db.load_full()
     assert len(rankings.index) == 5
     assert len(rankings.columns) == 22284
 
+
 def test_load(db, gs):
     rankings = db.load(gs)
     assert len(rankings.index) == 5
diff --git a/tests/test_featureseq.py b/tests/test_featureseq.py
index 0372eb9..28564e0 100644
--- a/tests/test_featureseq.py
+++ b/tests/test_featureseq.py
@@ -1,7 +1,7 @@
 # coding=utf-8
 
 
-from pyscenic.featureseq import Feature, FeatureSeq
+from pyscenic.featureseq import Feature
 
 
 def test_feature_overlap():
@@ -33,5 +33,3 @@ def test_feature_bp_overlap():
     assert f2.get_overlap_in_bp_with(f3) == 0
     assert f3.get_overlap_in_bp_with(f4) == 0
     assert f3.get_overlap_in_bp_with(f3) == len(f3)
-
-
diff --git a/tests/test_genesig.py b/tests/test_genesig.py
index 91c1b76..10c05b1 100644
--- a/tests/test_genesig.py
+++ b/tests/test_genesig.py
@@ -1,12 +1,10 @@
 # -*- coding: utf-8 -*-
 
-from pyscenic.genesig import GeneSignature, Regulon
-from configparser import ConfigParser
-import os
-import pytest
 import attr
+import pytest
 from pkg_resources import resource_filename
 
+from pyscenic.genesig import GeneSignature, Regulon
 
 TEST_SIGNATURE = "msigdb_cancer_c6"
 TEST_SIGNATURE_FNAME = resource_filename('resources.tests', "c6.all.v6.1.symbols.gmt")
@@ -40,7 +38,7 @@ def test_init3():
     assert len(gs1) == 2
     assert gs1.gene2weight['TP53'] == 0.5
     assert gs1.gene2weight['SOX4'] == 0.75
-    
+
 
 def test_immut():
     gs1 = GeneSignature(name="test1", gene2weight={'TP53': 0.5, 'SOX4': 0.75})
@@ -81,6 +79,7 @@ def test_union1():
     assert 'SOX2' in gsu
     assert len(gsu) == 3
 
+
 def test_union3():
     gs1 = GeneSignature(name="test1", gene2weight={'TP53': 0.8, 'SOX4': 0.75})
     gs2 = GeneSignature(name="test1", gene2weight={'TP53': 0.3, 'SOX2': 0.60})
@@ -134,8 +133,7 @@ def test_regulon():
 
 
 def test_load_gmt():
-    gss = GeneSignature.from_gmt(field_separator='\t', gene_separator='\t',
-                                 fname=TEST_SIGNATURE_FNAME)
+    gss = GeneSignature.from_gmt(field_separator='\t', gene_separator='\t', fname=TEST_SIGNATURE_FNAME)
     # http://software.broadinstitute.org/gsea/msigdb/collections.jsp#C6
     assert len(gss) == 189
     assert gss[0].name == "GLI1_UP.V1_DN"
@@ -144,8 +142,7 @@ def test_load_gmt():
 
 
 def test_add():
-    gss = GeneSignature.from_gmt(field_separator='\t', gene_separator='\t',
-                                 fname=TEST_SIGNATURE_FNAME)
+    gss = GeneSignature.from_gmt(field_separator='\t', gene_separator='\t', fname=TEST_SIGNATURE_FNAME)
     res = gss[0].add("MEF2")
     assert "MEF2" in res
 
@@ -162,6 +159,7 @@ def test_noweights():
     assert reg2['TP53'] == 1.0
     assert isinstance(reg2, Regulon)
 
+
 def test_head():
     gs1 = GeneSignature(name="test1", gene2weight={'TP53': 0.8, 'SOX4': 0.75})
     gs2 = gs1.head(1)
@@ -171,4 +169,3 @@ def test_head():
     assert gs2['TP53'] == 0.8
     assert gs2['SOX4'] == 0.75
     assert len(gs2) == 2
-
diff --git a/tests/test_math.py b/tests/test_math.py
index d6f3b4a..dea30c2 100644
--- a/tests/test_math.py
+++ b/tests/test_math.py
@@ -1,12 +1,13 @@
 # coding=utf-8
 
-import pytest
 import numpy as np
-from pyscenic.math import masked_rho, masked_rho_2d
+import pytest
 
+from pyscenic.math import masked_rho
 
-def test_masked_rho():
-    x = np.array([1, 2, -1, 0.0, 4], dtype=np.float)
-    y = np.array([1, 2, 0.0, -8, 4],  dtype=np.float)
-    assert pytest.approx(masked_rho(x, y, 0.0), 0.00001) == np.corrcoef(np.array([1.0, 2.0, 4.0]), np.array([1.0, 2.0, 4.0]))[0, 1]
 
+def test_masked_rho():
+    x = np.array([1, 2, -1, 0.0, 4], dtype=float)
+    y = np.array([1, 2, 0.0, -8, 4], dtype=float)
+    assert pytest.approx(masked_rho(x, y, 0.0), 0.00001) == \
+           np.corrcoef(np.array([1.0, 2.0, 4.0]), np.array([1.0, 2.0, 4.0]))[0, 1]
diff --git a/tests/test_parquetdb.py b/tests/test_parquetdb.py
new file mode 100644
index 0000000..21791f6
--- /dev/null
+++ b/tests/test_parquetdb.py
@@ -0,0 +1,50 @@
+# -*- coding: utf-8 -*-
+
+import pytest
+from pkg_resources import resource_filename
+
+from pyscenic.genesig import GeneSignature
+from pyscenic.rnkdb import ParquetRankingDatabase as RankingDatabase
+
+TEST_DATABASE_FNAME = resource_filename('resources.tests', "hg19-tss-centered-10kb-10species.mc9nr.parquet")
+TEST_DATABASE_NAME = "hg19-tss-centered-10kb-10species.mc9nr"
+TEST_SIGNATURE_FNAME = resource_filename('resources.tests', "c6.all.v6.1.symbols.gmt")
+
+##################################################
+# temporarily disable testing of the parqet databases until we have a working test database
+from os import path
+pytestmark = pytest.mark.skipif(not path.exists(TEST_DATABASE_FNAME), reason="Parquet testing is temporarily disabled.")
+##################################################
+
+@pytest.fixture
+def db():
+    return RankingDatabase(TEST_DATABASE_FNAME, TEST_DATABASE_NAME)
+
+
+@pytest.fixture
+def gs():
+    return GeneSignature.from_gmt(TEST_SIGNATURE_FNAME, gene_separator="\t", field_separator="\t", )[0]
+
+
+def test_init(db):
+    assert db.name == TEST_DATABASE_NAME
+
+
+def test_total_genes(db):
+    assert db.total_genes == 22284
+
+
+def test_genes(db):
+    assert len(db.genes) == 22284
+
+
+def test_load_full(db):
+    rankings = db.load_full()
+    assert len(rankings.index) == 5
+    assert len(rankings.columns) == 22284
+
+
+def test_load(db, gs):
+    rankings = db.load(gs)
+    assert len(rankings.index) == 5
+    assert len(rankings.columns) == 29
diff --git a/tests/test_recovery.py b/tests/test_recovery.py
index 322cb53..16fb7ae 100644
--- a/tests/test_recovery.py
+++ b/tests/test_recovery.py
@@ -1,13 +1,12 @@
 # -*- coding: utf-8 -*-
 
-from pyscenic.recovery import enrichment4features as enrichment, auc1d, weighted_auc1d, rcc2d
-
-import pytest
 import numpy as np
-from pyscenic.rnkdb import FeatherRankingDatabase as RankingDatabase
-from pyscenic.genesig import GeneSignature
+import pytest
 from pkg_resources import resource_filename
 
+from pyscenic.genesig import GeneSignature
+from pyscenic.recovery import enrichment4features as enrichment, auc1d, weighted_auc1d, rcc2d
+from pyscenic.rnkdb import FeatherRankingDatabase as RankingDatabase
 
 TEST_DATABASE_FNAME = resource_filename('resources.tests', "hg19-tss-centered-10kb-10species.mc9nr.feather")
 TEST_DATABASE_NAME = "hg19-tss-centered-10kb-10species.mc9nr"
@@ -16,13 +15,12 @@
 
 @pytest.fixture
 def db():
-    return RankingDatabase(TEST_DATABASE_FNAME, TEST_DATABASE_NAME,)
+    return RankingDatabase(TEST_DATABASE_FNAME, TEST_DATABASE_NAME, )
 
 
 @pytest.fixture
 def gs():
-    return GeneSignature.from_gmt(TEST_SIGNATURE_FNAME,
-                                  gene_separator="\t", field_separator="\t", )[0]
+    return GeneSignature.from_gmt(TEST_SIGNATURE_FNAME, gene_separator="\t", field_separator="\t", )[0]
 
 
 def test_enrichment(db, gs):
@@ -32,38 +30,48 @@ def test_enrichment(db, gs):
 def test_auc1d_1():
     # Check if AUC is calculated correctly when a gene is recovered at the rank threshold.
     # In the python implementation it should be included in the AUC calculation.
-    # For the current gene list the non-normalized AUC should be (2*1)+(2*2)+(2*3)+(1*4) = 16.
+    # For the current gene list the non-normalized AUC should be: (2*1) + (2*2) + (2*3) + (1*4) = 16.
+
     total_genes = 100
-    auc_rank_threshold = 8 # This rank threshold is inclusive.
-    ranking = np.asarray([2, 4, 6, 8]) - 1 # The databases have a zero-based ranking system
+    # This rank threshold is inclusive.
+    auc_rank_threshold = 8
+    # The databases have a zero-based ranking system.
+    ranking = np.asarray([2, 4, 6, 8]) - 1
     auc_max = auc_rank_threshold * total_genes
-    assert 16.0/800 == auc1d(ranking, auc_rank_threshold, auc_max)
+    assert 16.0 / 800 == auc1d(ranking, auc_rank_threshold, auc_max)
 
 
 def test_auc1d_2():
-    # For the current gene list the non-normalized AUC should be (2*1)+(2*2)+(3*3) = 15.
+    # For the current gene list the non-normalized AUC should be: (2*1) + (2*2) + (3*3) = 15.
+
     total_genes = 100
-    auc_rank_threshold = 8 # This rank threshold is inclusive.
-    ranking = np.asarray([2, 4, 6]) - 1 # The databases have a zero-based ranking system
+    # This rank threshold is inclusive.
+    auc_rank_threshold = 8
+    # The databases have a zero-based ranking system.
+    ranking = np.asarray([2, 4, 6]) - 1
     auc_max = auc_rank_threshold * total_genes
-    assert 15.0/800 == auc1d(ranking, auc_rank_threshold, auc_max)
+    assert 15.0 / 800 == auc1d(ranking, auc_rank_threshold, auc_max)
 
 
 def test_weighted_auc1d():
     # CAVE: In python the ranking databases are 0-based. The only difference a 1-based system has on the calc
-    # of the AUC is that in the latter the rank threshold would not be included. This has an infuence on the
+    # of the AUC is that in the latter the rank threshold would not be included. This has an influence on the
     # normalization factor max AUC but nothing else.
+
     total_genes = 100
     auc_rank_threshold = 8
-    ranking = np.asarray([2, 4, 6]) - 1 # The databases have a zero-based ranking system
+    # The databases have a zero-based ranking system.
+    ranking = np.asarray([2, 4, 6]) - 1
     weights = np.ones(len(ranking))
-    auc_max = 1.0 # Disable normalization.
+    # Disable normalization.
+    auc_max = 1.0
     assert 15.0 == weighted_auc1d(ranking, weights, auc_rank_threshold, auc_max)
 
 
 def test_weighted_auc1d_batch():
-    # The assumption taken here is that for weights uniformely being set to 1.0, auc1d and
-    # weighted_auc1d should always have the same output.
+    # The assumption taken here is that for weights uniformly being set to 1.0,
+    # auc1d and weighted_auc1d should always have the same output.
+
     N = 100000
     total_genes = 100
     for _ in range(N):
@@ -71,23 +79,28 @@ def test_weighted_auc1d_batch():
         n_genes = np.random.randint(20, high=total_genes)
         ranking = np.random.randint(low=0, high=total_genes, size=n_genes)
         weights = np.ones(n_genes)
-        auc_max = 1.0 # Disable normalization.
-        assert auc1d(ranking, auc_rank_threshold, auc_max) == weighted_auc1d(ranking, weights, auc_rank_threshold, auc_max)
+        # Disable normalization.
+        auc_max = 1.0
+        assert auc1d(ranking, auc_rank_threshold, auc_max) \
+               == weighted_auc1d(ranking, weights, auc_rank_threshold, auc_max)
 
 
 def test_weighted_rcc2d_batch():
-    # The assumption taken here is that for weights uniformely being set to 1.0, auc1d and
-    # weighted_auc1d should always have the same output.
+    # The assumption taken here is that for weights uniformly being set to 1.0,
+    # auc1d and weighted_auc1d should always have the same output.
+
     N = 100000
     total_genes = 100
     for _ in range(N):
-        #rccs[row_idx, :] = np.cumsum(np.bincount(curranking, weights=weights)[:rank_threshold])
-        #ValueError: could not broadcast input array from shape (89) into shape (97)
+        # rccs[row_idx, :] = np.cumsum(np.bincount(curranking, weights=weights)[:rank_threshold])
+        # ValueError: could not broadcast input array from shape (89) into shape (97)
         auc_rank_threshold = np.random.randint(2, high=total_genes)
         n_genes = np.random.randint(20, high=total_genes)
         ranking = np.random.randint(low=0, high=total_genes, size=n_genes)
         rankings = ranking.reshape((1, n_genes))
         rankings = np.append(rankings, np.full(shape=(1, 1), fill_value=total_genes), axis=1)
         weights = np.ones(n_genes)
-        auc_max = 1.0 # Disable normalization.
-        assert rcc2d(rankings, np.insert(weights, len(weights), 0.0), total_genes)[:, :auc_rank_threshold].sum(axis=1) == weighted_auc1d(ranking, weights, auc_rank_threshold, auc_max)
+        # Disable normalization.
+        auc_max = 1.0
+        assert rcc2d(rankings, np.insert(weights, len(weights), 0.0), total_genes)[:, :auc_rank_threshold].sum(axis=1) \
+               == weighted_auc1d(ranking, weights, auc_rank_threshold, auc_max)
diff --git a/tests/test_sqlitedb.py b/tests/test_sqlitedb.py
index 70b8431..ef6a9a6 100644
--- a/tests/test_sqlitedb.py
+++ b/tests/test_sqlitedb.py
@@ -1,11 +1,10 @@
 # -*- coding: utf-8 -*-
 
 import pytest
-from pyscenic.rnkdb import SQLiteRankingDatabase as RankingDatabase
-from pyscenic.genesig import GeneSignature
 from pkg_resources import resource_filename
 
-
+from pyscenic.genesig import GeneSignature
+from pyscenic.rnkdb import SQLiteRankingDatabase as RankingDatabase
 
 TEST_DATABASE_FNAME = resource_filename('resources.tests', "hg19-tss-centered-5kb-10species.mc9nr.db")
 TEST_DATABASE_NAME = "hg19-tss-centered-5kb-10species.mc9nr"
@@ -16,23 +15,30 @@
 def db():
     return RankingDatabase(TEST_DATABASE_FNAME, TEST_DATABASE_NAME)
 
+
 @pytest.fixture
 def gs():
-    return GeneSignature.from_gmt(TEST_SIGNATURE_FNAME,
-                                  gene_separator="\t", field_separator="\t", )[0]
+    return GeneSignature.from_gmt(TEST_SIGNATURE_FNAME, gene_separator="\t", field_separator="\t", )[0]
 
 
 def test_init(db):
     assert db.name == TEST_DATABASE_NAME
 
+
 def test_total_genes(db):
     assert db.total_genes == 29
 
+
+def test_genes(db):
+    assert len(db.genes) == 29
+
+
 def test_load_full(db):
     rankings = db.load_full()
     assert len(rankings.index) == 24453
     assert len(rankings.columns) == 29
 
+
 def test_load(db, gs):
     rankings = db.load(gs)
     assert len(rankings.index) == 24453