search

kevinrue / velociraptor

Toolkit for Single-Cell Velocity
https://kevinrue.github.io/velociraptor/
Other
53 stars 10 forks source link

Output of gridVectors contains NaN and NA values #26

Open jma1991 opened 3 years ago

jma1991 commented 3 years ago

Problem summary I am trying to annotate a UMAP plot with arrows corresponding to the RNA velocity vectors, as outlined in the RNA velocity section of the OSCA book. However, the gridVectors function returns a matrix with NaN or NA for all of the entries in the "end" columns. Unfortunately, I'm not sure if this is a problem with my data, scVelo itself, or velociraptor.

Reproducible example I have uploaded SingleCellExperiment objects of the expression data and the scvelo output (https://www.dropbox.com/sh/tx3lyej9925zo2w/AAAc5nqVJPY0gF9Paa_KGA7ra?dl=0). Below is some code to reproduce the problem:

> library(velociraptor)
> sce <- readRDS("SCE.rds")
> vel <- readRDS("scvelo.rds")
> mat <- embedVelocity(x = sce, vobj = vel, use.dimred = "UMAP")
> vec <- gridVectors(x = sce, embedded = mat, use.dimred = "UMAP")
> head(vec, n = 20)
   start.UMAP.1 start.UMAP.2 end.UMAP.1 end.UMAP.2
1   -2.87237923 -1.664514254        NaN        NaN
2   -0.52719200 -3.804835333        NaN        NaN
3    1.82489239 -6.936517152        NaN        NaN
4    2.44685568 -4.877971456        NaN        NaN
5   -2.82688955 -3.771670043        NaN        NaN
6    2.75407094 -6.446214525        NaN        NaN
7   -0.49538154 -1.080135230        NaN        NaN
8   -1.15055107 -5.945217601        NaN        NaN
10   0.47848420 -1.764962397        NaN        NaN
11  -2.81111142 -3.279491583        NaN        NaN
12  -1.05147457 -4.382347539        NaN        NaN
13   1.76235273 -0.002760358        NaN        NaN
14 -15.20460842  2.421717005        NaN        NaN
15   1.39724138 -1.250904078        NaN        NaN
16  -3.49125886 -3.246118226        NaN        NaN
17  -0.25546962  2.671420269         NA         NA
18  -0.61634884 -4.323074135        NaN        NaN
19   4.01094177 -2.191837065        NaN        NaN
20   0.07210263  2.553780907         NA         NA
21   0.08831714 -1.710850282        NaN        NaN

However, when I set the scale parameter to FALSE the NA values are replaced by numeric values:

> vec <- gridVectors(x = sce, embedded = mat, use.dimred = "UMAP", scale = FALSE)
> head(vec, n = 20)
   start.UMAP.1 start.UMAP.2  end.UMAP.1 end.UMAP.2
1   -2.87237923 -1.664514254         NaN        NaN
2   -0.52719200 -3.804835333         NaN        NaN
3    1.82489239 -6.936517152         NaN        NaN
4    2.44685568 -4.877971456         NaN        NaN
5   -2.82688955 -3.771670043         NaN        NaN
6    2.75407094 -6.446214525         NaN        NaN
7   -0.49538154 -1.080135230         NaN        NaN
8   -1.15055107 -5.945217601         NaN        NaN
10   0.47848420 -1.764962397         NaN        NaN
11  -2.81111142 -3.279491583         NaN        NaN
12  -1.05147457 -4.382347539         NaN        NaN
13   1.76235273 -0.002760358         NaN        NaN
14 -15.20460842  2.421717005         NaN        NaN
15   1.39724138 -1.250904078         NaN        NaN
16  -3.49125886 -3.246118226         NaN        NaN
17  -0.25546962  2.671420269 -0.19771867   2.237185
18  -0.61634884 -4.323074135         NaN        NaN
19   4.01094177 -2.191837065         NaN        NaN
20   0.07210263  2.553780907  0.08957561   2.221970
21   0.08831714 -1.710850282         NaN        NaN

Session information

R version 4.0.3 (2020-10-10)
Platform: x86_64-apple-darwin13.4.0 (64-bit)
Running under: macOS Catalina 10.15.6

Matrix products: default
BLAS/LAPACK: /Users/james/opt/miniconda3/envs/orchestra/lib/libopenblasp-r0.3.12.dylib

locale:
[1] en_GB.UTF-8/en_GB.UTF-8/en_GB.UTF-8/C/en_GB.UTF-8/en_GB.UTF-8

attached base packages:
[1] parallel  stats4    stats     graphics  grDevices utils     datasets
[8] methods   base

other attached packages:
 [1] velociraptor_1.0.0          SummarizedExperiment_1.20.0
 [3] Biobase_2.50.0              GenomicRanges_1.42.0
 [5] GenomeInfoDb_1.26.0         IRanges_2.24.0
 [7] S4Vectors_0.28.0            BiocGenerics_0.36.0
 [9] MatrixGenerics_1.2.0        matrixStats_0.57.0

loaded via a namespace (and not attached):
 [1] Rcpp_1.0.5                  zellkonverter_1.0.0
 [3] BiocSingular_1.6.0          XVector_0.30.0
 [5] zlibbioc_1.36.0             rappdirs_0.3.1
 [7] BiocParallel_1.24.0         lattice_0.20-41
 [9] sparseMatrixStats_1.2.0     DelayedMatrixStats_1.12.0
[11] tools_4.0.3                 rsvd_1.0.3
[13] grid_4.0.3                  scuttle_1.0.0
[15] irlba_2.3.3                 basilisk_1.2.0
[17] Matrix_1.2-18               GenomeInfoDbData_1.2.4
[19] bitops_1.0-6                basilisk.utils_1.2.0
[21] SingleCellExperiment_1.12.0 RCurl_1.98-1.2
[23] DelayedArray_0.16.0         compiler_4.0.3
[25] filelock_1.0.2              beachmat_2.6.0
[27] jsonlite_1.7.1              reticulate_1.18
kevinrue commented 3 years ago

Hi @jma1991

I've run your code with the data from Dropbox.

From what I can see, mat already contains a large number of NaN, so the issues seems to lie in embedVelocity rather than gridVectors.

embedVelocity is not doing much more than running scvelo under the hood, so my guess is that the NaN are generated by scvelo rather than velociraptor.

I would suggest running your data through scvelo in a Python session, to confirm that is the case.

As a side note, gridVectors(..., scale = FALSE) is just avoiding an additional computation whereby the NaN interfere with the scaling of the remaining valid numeric values, turning those into NAs. https://github.com/kevinrue/velociraptor/blob/49fbf28b175194f07252e86d3a607f4915872dd6/R/gridVectors.R#L69 As I pointed out above, this is not the source of the issue, but rather a consequence of the NaN introduced earlier during the embedVelocity call.

jma1991 commented 3 years ago

Thanks for the quick reply, I will have a go at using scVelo from within Python and report back.

Hi @jma1991

I've run your code with the data from Dropbox.

From what I can see, mat already contains a large number of NaN, so the issues seems to lie in embedVelocity rather than gridVectors.

embedVelocity is not doing much more than running scvelo under the hood, so my guess is that the NaN are generated by scvelo rather than velociraptor.

I would suggest running your data through scvelo in a Python session, to confirm that is the case.

As a side note, gridVectors(..., scale = FALSE) is just avoiding an additional computation whereby the NaN interfere with the scaling of the remaining valid numeric values, turning those into NAs.

https://github.com/kevinrue/velociraptor/blob/49fbf28b175194f07252e86d3a607f4915872dd6/R/gridVectors.R#L69

As I pointed out above, this is not the source of the issue, but rather a consequence of the NaN introduced earlier during the embedVelocity call.

jma1991 commented 3 years ago

Hi @kevinrue

I've managed to run scVelo from within Python and everything seems to be working fine (i.e. it generates RNA velocity vectors for most cells and I can plot them on the UMAP embedding, see linked image: https://ibb.co/zJs4t04). To run this test I first activated the velociraptor environment which gets installed when I first run the scvelo function. From within this environment I then imported the SCE object, converted it to an AnnData object and ran the usual scvelo commands (I tried to match the commands velociraptor sends to scvelo internally):

import scvelo as scv

import anndata2ri

from rpy2.robjects import r

anndata2ri.activate()

scv.set_figure_params()

r('library(SingleCellExperiment)')

adata = r('readRDS("SCE.rds")')

hvg = adata.var['HVG'] > 0 # turn highly variable gene series back into a boolean

sdata = adata[:, hvg]

scv.pp.moments(sdata)

scv.tl.velocity(sdata, mode = "steady_state")

scv.tl.velocity_graph(sdata)

scv.tl.velocity_pseudotime(sdata)

scv.tl.velocity_confidence(sdata)

scv.pl.velocity_embedding(sdata, basis = 'X_umap')

I've uploaded the original (adata.obj) and subsetted (sdata.obj) AnnData objects to the same Dropbox folder in case you need to take a look at them yourself. Below is a list of packages in the conda environment:

$ conda list
# packages in environment at /Users/james/Library/Caches/basilisk/1.2.0/velociraptor-1.0.0/env:
#
# Name                    Version                   Build  Channel
alabaster                 0.7.12                     py_0    conda-forge
anndata                   0.7.4            py37hf985489_1    conda-forge
babel                     2.8.0                      py_0    conda-forge
blosc                     1.20.1               hb1e8313_0    conda-forge
brotlipy                  0.7.0           py37h395d20d_1001    conda-forge
bzip2                     1.0.8                haf1e3a3_3    conda-forge
c-ares                    1.11.0                        0    bioconda
ca-certificates           2020.11.8            h033912b_0    conda-forge
certifi                   2020.6.20        py37h2987424_2    conda-forge
cffi                      1.14.1           py37hed5b41f_0
chardet                   3.0.4           py37h2987424_1008    conda-forge
colorama                  0.4.3                      py_0    conda-forge
cryptography              3.1.1            py37h55b06a8_1    conda-forge
cycler                    0.10.0                     py_2    conda-forge
decorator                 4.4.2                      py_0    conda-forge
docutils                  0.16             py37h2987424_2    conda-forge
freetype                  2.10.4               h3f75d11_0    conda-forge
get_version               2.1                        py_1    conda-forge
h5py                      2.10.0          nompi_py37h6d80030_105    conda-forge
hdf5                      1.10.6          nompi_h2ccf146_1110    conda-forge
idna                      2.10               pyh9f0ad1d_0    conda-forge
imagesize                 1.2.0                      py_0    conda-forge
importlib-metadata        1.7.0            py37hc8dfbb8_0    conda-forge
importlib_metadata        1.7.0                         0    conda-forge
jinja2                    2.11.2             pyh9f0ad1d_0    conda-forge
joblib                    0.16.0                     py_0    conda-forge
jpeg                      9d                   hbcb3906_0    conda-forge
kiwisolver                1.2.0            py37hb449ec0_1    conda-forge
krb5                      1.17.1               h75d18d8_3    conda-forge
lcms2                     2.11                 h11f7e16_1    conda-forge
legacy-api-wrap           1.2                        py_0    conda-forge
libblas                   3.9.0                2_openblas    conda-forge
libcblas                  3.9.0                2_openblas    conda-forge
libcurl                   7.71.1               h9bf37e3_8    conda-forge
libcxx                    11.0.0               h439d374_0    conda-forge
libedit                   3.1.20191231         h0678c8f_2    conda-forge
libev                     4.33                 haf1e3a3_1    conda-forge
libffi                    3.3                  hb1e8313_2
libgfortran               5.0.0               h7cc5361_13    conda-forge
libgfortran5              9.3.0               h7cc5361_13    conda-forge
liblapack                 3.9.0                2_openblas    conda-forge
libllvm10                 10.0.1               h009f743_3    conda-forge
libnghttp2                1.41.0               h8a08a2b_1    conda-forge
libopenblas               0.3.12          openmp_h54245bb_1    conda-forge
libpng                    1.6.37               h7cec526_2    conda-forge
libssh2                   1.9.0                h8a08a2b_5    conda-forge
libtiff                   4.1.0                hca7d577_6    conda-forge
libwebp-base              1.1.0                hbcb3906_3    conda-forge
llvm-openmp               11.0.0               h73239a0_1    conda-forge
llvmlite                  0.34.0           py37h3986384_2    conda-forge
loompy                    2.0.16                     py_0    bioconda
lz4-c                     1.9.2                hb1e8313_3    conda-forge
markupsafe                1.1.1            py37h395d20d_2    conda-forge
matplotlib                3.3.1                         1    conda-forge
matplotlib-base           3.3.1            py37h886f89f_1    conda-forge
mock                      4.0.2            py37hc8dfbb8_1    conda-forge
natsort                   7.0.1                      py_0    conda-forge
ncurses                   6.2                  h2e338ed_3    conda-forge
networkx                  2.5                        py_0    conda-forge
numba                     0.51.2           py37h6d0141a_0    conda-forge
numexpr                   2.7.1            py37h6d0141a_3    conda-forge
numpy                     1.19.1           py37h1efc2f6_2    conda-forge
olefile                   0.46               pyh9f0ad1d_1    conda-forge
openssl                   1.1.1h               haf1e3a3_0    conda-forge
packaging                 20.4               pyh9f0ad1d_0    conda-forge
pandas                    1.1.2            py37hdadc0f0_0    conda-forge
patsy                     0.5.1                      py_0    conda-forge
pillow                    7.2.0            py37hf860fee_2    conda-forge
pip                       20.2.3                     py_0    conda-forge
pycparser                 2.20               pyh9f0ad1d_2    conda-forge
pygments                  2.7.0                      py_0    conda-forge
pyopenssl                 19.1.0                     py_1    conda-forge
pyparsing                 2.4.7              pyh9f0ad1d_0    conda-forge
pysocks                   1.7.1            py37h2987424_2    conda-forge
pytables                  3.6.1            py37he95298d_3    conda-forge
python                    3.7.7                hf48f09d_4
python-dateutil           2.8.1                      py_0    conda-forge
python_abi                3.7                     1_cp37m    conda-forge
pytz                      2020.1             pyh9f0ad1d_0    conda-forge
readline                  8.0                  h0678c8f_2    conda-forge
requests                  2.24.0             pyh9f0ad1d_0    conda-forge
scanpy                    1.6.0                      py_0    bioconda
scikit-learn              0.23.2           py37hd891755_2    conda-forge
scipy                     1.5.2            py37h821cff1_2    conda-forge
scvelo                    0.2.2                      py_1    bioconda
seaborn                   0.11.0               h694c41f_1    conda-forge
seaborn-base              0.11.0             pyhd8ed1ab_1    conda-forge
setuptools                49.6.0           py37h2987424_2    conda-forge
setuptools-scm            4.1.2              pyh9f0ad1d_0    conda-forge
setuptools_scm            4.1.2                         0    conda-forge
sinfo                     0.3.1                      py_0    conda-forge
six                       1.15.0             pyh9f0ad1d_0    conda-forge
snowballstemmer           2.0.0                      py_0    conda-forge
sphinx                    3.2.1                      py_0    conda-forge
sphinxcontrib-applehelp   1.0.2                      py_0    conda-forge
sphinxcontrib-devhelp     1.0.2                      py_0    conda-forge
sphinxcontrib-htmlhelp    1.0.3                      py_0    conda-forge
sphinxcontrib-jsmath      1.0.1                      py_0    conda-forge
sphinxcontrib-qthelp      1.0.3                      py_0    conda-forge
sphinxcontrib-serializinghtml 1.1.4                      py_0    conda-forge
sqlite                    3.33.0               h960bd1c_1    conda-forge
statsmodels               0.12.0           py37h57c32b8_1    conda-forge
stdlib-list               0.6.0                    py37_0    conda-forge
tbb                       2019.9               ha1b3eb9_1    conda-forge
threadpoolctl             2.1.0              pyh5ca1d4c_0    conda-forge
tk                        8.6.10               hb0a8c7a_1    conda-forge
toml                      0.10.1             pyh9f0ad1d_0    conda-forge
tornado                   6.0.4            py37h60d8a13_2    conda-forge
tqdm                      4.49.0             pyh9f0ad1d_0    conda-forge
typing                    3.7.4.3          py37hc8dfbb8_1    conda-forge
umap-learn                0.4.6            py37hc8dfbb8_0    conda-forge
urllib3                   1.25.10                    py_0    conda-forge
wheel                     0.35.1             pyh9f0ad1d_0    conda-forge
xz                        5.2.5                haf1e3a3_1    conda-forge
zipp                      3.1.0                      py_0    conda-forge
zlib                      1.2.11            h7795811_1010    conda-forge
zstd                      1.4.5                h289c70a_2    conda-forge
(/Users/james/Library/Caches/basilisk/1.2.0/velociraptor-1.0.0/env)
LTLA commented 3 years ago

I can confirm that setting use.theirs=TRUE in scvelo() gets rid of the NA's. It also eats up all the RAM on my machine, but it does eventually return non-NA values in mat.

When I was running scvelo(sce, subset.row=rowData(sce)$HVG), I noticed:

computing velocities
WARNING: Too few genes are selected as velocity genes. Consider setting a lower threshold for min_r2 or min_likelihood.

It's a very chatty function - we haven't silenced the Python output, so you get all the noise in its full glory. This makes it hard to spot the WARNINGS when they do actually show up. But that warning does sound like it could be problematic.

Well, okay, I just bumped up the number of genes. I don't know how you get HVG, but this is what I did:

dec <- scran::modelGeneVar(sce)
vel <- scvelo(sce, subset.row=getTopHVGs(dec, n=5000))

I didn't check that the trend was sensible or anything, so you might not want to do exactly what I did. But it doesn't really matter, the point is that you can just ask for more genes from whatever algorithm that you're using to detect HVGs. Then hopefully the NAs will go away - they did for me, when I checked mat and vec.

Why? I don't really know, but if I had to guess, genes that are informative for characterizing variation are not necessarily informative for computing velocities, e.g., if a HVG doesn't have any intronic counts then it's pretty useless. I'd speculate that some cells didn't have enough counts for the few informative genes and just ended up with an unknown velocity.

jma1991 commented 3 years ago

Following your advice Aaron, I have loosened my threshold for selecting highly variable genes and have been able to generate RNA velocity vectors for nearly all of the cells. Hurray!

I'm still not sure why the previous number of highly variable genes "worked" within scVelo but not within velociraptor. I think I used the same commands for both runs, also skipping the "filter_and_normalize" option within scVelo. They should be the same?

I can confirm that setting use.theirs=TRUE in scvelo() gets rid of the NA's. It also eats up all the RAM on my machine, but it does eventually return non-NA values in mat.

When I was running scvelo(sce, subset.row=rowData(sce)$HVG), I noticed:

computing velocities
WARNING: Too few genes are selected as velocity genes. Consider setting a lower threshold for min_r2 or min_likelihood.

It's a very chatty function - we haven't silenced the Python output, so you get all the noise in its full glory. This makes it hard to spot the WARNINGS when they do actually show up. But that warning does sound like it could be problematic.

Well, okay, I just bumped up the number of genes. I don't know how you get HVG, but this is what I did:

dec <- scran::modelGeneVar(sce)
vel <- scvelo(sce, subset.row=getTopHVGs(dec, n=5000))

I didn't check that the trend was sensible or anything, so you might not want to do exactly what I did. But it doesn't really matter, the point is that you can just ask for more genes from whatever algorithm that you're using to detect HVGs. Then hopefully the NAs will go away - they did for me, when I checked mat and vec.

Why? I don't really know, but if I had to guess, genes that are informative for characterizing variation are not necessarily informative for computing velocities, e.g., if a HVG doesn't have any intronic counts then it's pretty useless. I'd speculate that some cells didn't have enough counts for the few informative genes and just ended up with an unknown velocity.

sunshine1126 commented 3 years ago

@kevinrue hello, I also encountered a similar problem there was many NA values when I used gridVectors. The data was processed refer to the website (https://jef.works/blog/2020/08/25/using-scvelo-in-R-using-reticulate/). My code is as follows：

library(reticulate)
library(tidyverse)
library(anndata)
library(velociraptor)
os <- import("os")
numpy <- import('numpy')
scv <- import("scvelo")
## scvelo
adata <- scv$datasets$pancreas()
adata

scv$pp$filter_and_normalize(adata) ## filter and normalize
scv$pp$moments(adata) ## compute moments
scv$tl$velocity(adata, mode='dynamical')
scv$tl$velocity_graph(adata)
scv$pl$velocity_embedding(adata, basis='umap',arrow_length=3, arrow_size=2, dpi=120) ## add velocity_umap

dat1 <- cbind(adata$obsm$X_umap,adata$obsm$velocity_umap)
dat1 <- as.data.frame(dat1)
out1 <- gridVectors(dat1[,1:2],dat1[,3:4],scale=F)

dat1

image

out1

image I'm not sure if this is a problem with my data, scVelo itself, or velociraptor. Can you help me solve it? Thanks！

csoneson commented 3 years ago

👋🏻 @sunshine1126 - I ran the same code as you provide above, and I am not able to reproduce the NA values; I get 

> head(out1)
    start.V1    start.V2     end.V1     end.V2
1   5.977051 -0.05535176   5.969925 -0.0175779
2  -9.968773  0.19162985  -9.979441  0.1761117
3   7.660351  0.48652162   7.658802  0.4883735
4 -11.144027  4.15642950 -11.122134  4.1894750
5   1.523521 -4.62369369   1.613057 -4.6157716
6 -11.033881  0.80893791 -11.030386  0.7857211

Could you check which version of the software packages you are using? E.g., the current release version of velociraptor is v1.2.0 (which bundles scVelo v0.2.2).

sunshine1126 commented 3 years ago

@csoneson Thanks for your reply! I check the current software version is as follows：

R version 4.0.2 (2020-06-22)
Platform: x86_64-pc-linux-gnu (64-bit)
Running under: Debian GNU/Linux 10 (buster)

other attached packages:
 [1] plotly_4.9.2.1              scater_1.16.2               Matrix_1.2-18              
 [4] ggthemes_4.2.0              ggpubr_0.2.5                magrittr_2.0.1             
 [7] REdaS_0.9.3                 circular_0.4-93             scran_1.16.0               
[10] basilisk.utils_1.2.2        basilisk_1.2.1              scuttle_1.0.4              
[13] celldex_1.0.0               scRNAseq_2.4.0              SingleCellExperiment_1.12.0
[16] SingleR_1.4.1               future_1.17.0               vroom_1.4.0                
[19] configr_0.3.5               futile.logger_1.4.3         shiny_1.4.0.2              
[22] BiocManager_1.30.15         corrplot_0.88               networkD3_0.4              
[25] viridis_0.5.1               viridisLite_0.3.0           ggstance_0.3.4             
[28] ggtree_2.4.2                forcats_0.5.0               stringr_1.4.0              
[31] purrr_0.3.4                 readr_1.3.1                 tidyr_1.0.2                
[34] tibble_3.0.1                tidyverse_1.3.0             phyloseq_1.34.0            
[37] ggdendro_0.1-20             glue_1.4.2                  openxlsx_4.1.4             
[40] patchwork_1.0.0             reshape2_1.4.4              ggplot2_3.3.0              
[43] dplyr_0.8.5                 plyr_1.8.6                  tensorflow_2.0.0           
[46] reticulate_1.15             velociraptor_1.0.0          SummarizedExperiment_1.20.0
[49] Biobase_2.50.0              GenomicRanges_1.42.0        GenomeInfoDb_1.26.7        
[52] IRanges_2.24.1              S4Vectors_0.28.1            BiocGenerics_0.36.1        
[55] MatrixGenerics_1.2.1        matrixStats_0.59.0          devtools_2.3.0             
[58] usethis_1.6.0              

loaded via a namespace (and not attached):
  [1] rappdirs_0.3.1                rtracklayer_1.50.0           
  [3] bit64_4.0.5                   irlba_2.3.3                  
  [5] DelayedArray_0.16.3           data.table_1.14.0            
  [7] RCurl_1.98-1.3                AnnotationFilter_1.14.0      
  [9] generics_0.1.0                GenomicFeatures_1.42.3       
 [11] callr_3.4.3                   lambda.r_1.2.4               
 [13] RSQLite_2.2.0                 bit_4.0.4                    
 [15] xml2_1.3.2                    lubridate_1.7.8              
 [17] httpuv_1.5.2                  assertthat_0.2.1             
 [19] hms_0.5.3                     promises_1.1.0               
 [21] fansi_0.5.0                   progress_1.2.2               
 [23] dbplyr_1.4.3                  readxl_1.3.1                 
 [25] igraph_1.2.5                  DBI_1.1.1                    
 [27] htmlwidgets_1.5.1             ellipsis_0.3.2               
 [29] backports_1.2.1               permute_0.9-5                
 [31] biomaRt_2.46.3                sparseMatrixStats_1.2.1      
 [33] vctrs_0.2.4                   remotes_2.1.1                
 [35] ensembldb_2.14.1              withr_2.2.0                  
 [37] vegan_2.5-7                   GenomicAlignments_1.26.0     
 [39] treeio_1.14.4                 prettyunits_1.1.1            
 [41] zellkonverter_1.0.3           cluster_2.1.0                
 [43] ExperimentHub_1.16.1          ape_5.5                      
 [45] lazyeval_0.2.2                crayon_1.4.1                 
 [47] edgeR_3.30.3                  pkgconfig_2.0.3              
 [49] vipor_0.4.5                   nlme_3.1-148                 
 [51] pkgload_1.0.2                 ProtGenerics_1.22.0          
 [53] rlang_0.4.5                   globals_0.14.0               
 [55] lifecycle_0.2.0               filelock_1.0.2               
 [57] BiocFileCache_1.14.0          modelr_0.1.6                 
 [59] rsvd_1.0.3                    AnnotationHub_2.22.1         
 [61] cellranger_1.1.0              rprojroot_1.3-2              
 [63] aplot_0.0.6                   boot_1.3-25                  
 [65] Rhdf5lib_1.12.1               reprex_0.3.0                 
 [67] RcppTOML_0.1.7                base64enc_0.1-3              
 [69] beeswarm_0.2.3                whisker_0.4                  
 [71] processx_3.4.2                ini_0.3.1                    
 [73] bitops_1.0-7                  rhdf5filters_1.2.1           
 [75] Biostrings_2.58.0             blob_1.2.1                   
 [77] DelayedMatrixStats_1.12.3     ggsignif_0.6.0               
 [79] beachmat_2.6.4                scales_1.1.0                 
 [81] memoise_1.1.0                 zlibbioc_1.36.0              
 [83] compiler_4.0.2                dqrng_0.2.1                  
 [85] Rsamtools_2.6.0               cli_2.0.2                    
 [87] ade4_1.7-16                   XVector_0.30.0               
 [89] listenv_0.8.0                 ps_1.6.0                     
 [91] formatR_1.11                  MASS_7.3-51.6                
 [93] mgcv_1.8-31                   tidyselect_1.0.0             
 [95] stringi_1.4.6                 yaml_2.2.1                   
 [97] BiocSingular_1.6.0            askpass_1.1                  
 [99] locfit_1.5-9.4                tools_4.0.2                  
[101] rstudioapi_0.11               foreach_1.5.1                
[103] gridExtra_2.3                 digest_0.6.27                
[105] rvcheck_0.1.8                 Rcpp_1.0.6                   
[107] broom_0.5.6                   BiocVersion_3.12.0           
[109] later_1.2.0                   httr_1.4.1                   
[111] AnnotationDbi_1.52.0          colorspace_2.0-1             
[113] rvest_0.3.5                   XML_3.99-0.3                 
[115] fs_1.5.0                      splines_4.0.2                
[117] statmod_1.4.34                tidytree_0.3.4               
[119] multtest_2.46.0               sessioninfo_1.1.1            
[121] xtable_1.8-4                  jsonlite_1.7.2               
[123] futile.options_1.0.1          testthat_2.3.2               
[125] R6_2.4.1                      pillar_1.4.3                 
[127] htmltools_0.4.0               mime_0.10                    
[129] fastmap_1.1.0                 BiocParallel_1.24.1          
[131] BiocNeighbors_1.8.2           interactiveDisplayBase_1.28.0
[133] codetools_0.2-16              mvtnorm_1.1-0                
[135] pkgbuild_1.0.6                lattice_0.20-41              
[137] curl_4.3.1                    ggbeeswarm_0.6.0             
[139] tfruns_1.4                    zip_2.0.4                    
[141] openssl_1.4.4                 survival_3.1-12              
[143] limma_3.44.3                  desc_1.2.0                   
[145] biomformat_1.18.0             munsell_0.5.0                
[147] rhdf5_2.34.0                  GenomeInfoDbData_1.2.4       
[149] iterators_1.0.13              haven_2.2.0                  
[151] gtable_0.3.0                 
# packages in environment at ~/miniconda3/envs/r-velocity:
#
# Name                    Version                   Build  Channel
_libgcc_mutex             0.1                 conda_forge    conda-forge
_openmp_mutex             4.5                       1_gnu    conda-forge
anndata                   0.7.6            py38h578d9bd_0    conda-forge
bzip2                     1.0.8                h7f98852_4    conda-forge
c-ares                    1.17.1               h7f98852_1    conda-forge
ca-certificates           2021.5.30            ha878542_0    conda-forge
cached-property           1.5.2                hd8ed1ab_1    conda-forge
cached_property           1.5.2              pyha770c72_1    conda-forge
certifi                   2021.5.30        py38h578d9bd_0    conda-forge
click                     8.0.1                    pypi_0    pypi
cycler                    0.10.0                   pypi_0    pypi
decorator                 4.4.2                    pypi_0    pypi
dunamai                   1.5.5                    pypi_0    pypi
get-version               3.0                      pypi_0    pypi
h5py                      3.2.1           nompi_py38h9915d05_100    conda-forge
hdf5                      1.10.6          nompi_h6a2412b_1114    conda-forge
htslib                    1.12                 h9093b5e_1    bioconda
joblib                    1.0.1                    pypi_0    pypi
kiwisolver                1.3.1                    pypi_0    pypi
krb5                      1.19.1               hcc1bbae_0    conda-forge
ld_impl_linux-64          2.35.1               hea4e1c9_2    conda-forge
legacy-api-wrap           1.2                      pypi_0    pypi
libblas                   3.9.0                9_openblas    conda-forge
libcblas                  3.9.0                9_openblas    conda-forge
libcurl                   7.77.0               h2574ce0_0    conda-forge
libdeflate                1.7                  h7f98852_5    conda-forge
libedit                   3.1.20191231         he28a2e2_2    conda-forge
libev                     4.33                 h516909a_1    conda-forge
libffi                    3.3                  h58526e2_2    conda-forge
libgcc-ng                 9.3.0               h2828fa1_19    conda-forge
libgfortran-ng            9.3.0               hff62375_19    conda-forge
libgfortran5              9.3.0               hff62375_19    conda-forge
libgomp                   9.3.0               h2828fa1_19    conda-forge
liblapack                 3.9.0                9_openblas    conda-forge
libnghttp2                1.43.0               h812cca2_0    conda-forge
libopenblas               0.3.15          pthreads_h8fe5266_1    conda-forge
libssh2                   1.9.0                ha56f1ee_6    conda-forge
libstdcxx-ng              9.3.0               h6de172a_19    conda-forge
llvmlite                  0.36.0                   pypi_0    pypi
loompy                    3.0.6                    pypi_0    pypi
matplotlib                3.4.2                    pypi_0    pypi
natsort                   7.1.1              pyhd8ed1ab_0    conda-forge
ncurses                   6.2                  h58526e2_4    conda-forge
networkx                  2.5.1                    pypi_0    pypi
numba                     0.53.1                   pypi_0    pypi
numexpr                   2.7.3                    pypi_0    pypi
numpy                     1.21.0           py38h9894fe3_0    conda-forge
numpy-groupies            0.9.13                   pypi_0    pypi
openssl                   1.1.1k               h7f98852_0    conda-forge
packaging                 20.9               pyh44b312d_0    conda-forge
pandas                    1.2.4                    pypi_0    pypi
patsy                     0.5.1                    pypi_0    pypi
pillow                    8.2.0                    pypi_0    pypi
pip                       21.1.2             pyhd8ed1ab_0    conda-forge
pyparsing                 2.4.7              pyh9f0ad1d_0    conda-forge
python                    3.8.10          h49503c6_1_cpython    conda-forge
python-dateutil           2.8.1                      py_0    conda-forge
python_abi                3.8                      1_cp38    conda-forge
pytz                      2021.1             pyhd8ed1ab_0    conda-forge
readline                  8.1                  h46c0cb4_0    conda-forge
samtools                  1.12                 h9aed4be_1    bioconda
scanpy                    1.7.2                    pypi_0    pypi
scikit-learn              0.24.2                   pypi_0    pypi
scipy                     1.6.3            py38h7b17777_0    conda-forge
scvelo                    0.2.3                    pypi_0    pypi
seaborn                   0.11.1                   pypi_0    pypi
setuptools                49.6.0           py38h578d9bd_3    conda-forge
sinfo                     0.3.4                    pypi_0    pypi
six                       1.16.0             pyh6c4a22f_0    conda-forge
sqlite                    3.36.0               h9cd32fc_0    conda-forge
statsmodels               0.12.2                   pypi_0    pypi
stdlib-list               0.8.0                    pypi_0    pypi
tables                    3.6.1                    pypi_0    pypi
threadpoolctl             2.1.0                    pypi_0    pypi
tk                        8.6.10               h21135ba_1    conda-forge
tqdm                      4.61.0                   pypi_0    pypi
umap-learn                0.4.6                    pypi_0    pypi
wheel                     0.36.2             pyhd3deb0d_0    conda-forge
xlrd                      1.2.0                    pypi_0    pypi
xz                        5.2.5                h516909a_1    conda-forge
zlib                      1.2.11            h516909a_1010    conda-forge

I found that our current release version is not the same，including velociraptor and scvelo.

csoneson commented 3 years ago

I would suggest to first update velociraptor to the current version. Then, check whether there are any NAs already in dat1 (in which case it comes from scVelo, since velociraptor was not used to get to this point - maybe too few genes are being selected like in the previous issue above). If there are no NAs in dat1, I think you can step through the gridVectors() function step by step and see where the issue occurs (code is here).