-
Notifications
You must be signed in to change notification settings - Fork 10
New issue
Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.
By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.
Already on GitHub? Sign in to your account
Output of gridVectors contains NaN and NA values #26
Comments
Hi @jma1991 I've run your code with the data from Dropbox. From what I can see,
I would suggest running your data through As a side note, Line 69 in 49fbf28
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.
|
Thanks for the quick reply, I will have a go at using scVelo from within Python and report back.
|
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) |
I can confirm that setting When I was running
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 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 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. |
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?
|
@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/).
|
👋🏻 @sunshine1126 - I ran the same code as you provide above, and I am not able to reproduce the NA values; I get
Could you check which version of the software packages you are using? E.g., the current release version of |
@csoneson Thanks for your reply! I check the current software version is as follows:
I found that our current release version is not the same,including velociraptor and scvelo. |
I would suggest to first update |
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:
However, when I set the scale parameter to FALSE the NA values are replaced by numeric values:
Session information
The text was updated successfully, but these errors were encountered: