Comprehensive RNA velocity by modeling the cascade of gene regulation, transcription, and splicing from single-cell RNA sequencing data with TSvelo
- State Key Laboratory of Biopharmaceutical Preparation and Delivery, Institute of Process Engineering, Chinese Academy of Sciences, China
- Institute of Image Processing and Pattern Recognition, Shanghai Jiao Tong University, and Key Laboratory of System Control and Information Processing, Ministry of Education of China, China
Figures
Figure 1 with 7 supplements
The framework of TSvelo.
(a) The preprocessing strategy in TSvelo, including velocity genes selection and initial state detection. (b) The neural ordinary differential equation (ODE) model and its optimization in TSvelo, where the parameters and latent time are optimized iteratively. The downstream application tasks of TSvelo, including precise transcription–unspliced–spliced 3D phase portrait fitting (c), cell fate prediction using predicted RNA velocity (d), and gene expression pattern analysis for multi-lineage dataset (e).
Figure 1—figure supplement 1
Initial state detection and the 2D velocity stream obtained by TSvelo on multiple simulation datasets.
Figure 1—figure supplement 2
Comparison between TSvelo and baseline approaches using the simulation data.
For these boxplots in the top row, the whiskers extend to 1.5× the interquartile range of the distribution from the hinge. The sample sizes are 600, 1200 and 1800 for one branch, tow branches and three branches respectively.
Figure 1—figure supplement 3
The U-to-S delay scores for each Leiden cluster when treated as the initial state across all datasets.
Figure 1—figure supplement 4
The quantitative comparison of TSvelo modeling results on the pancreas dataset using different transcription factor (TF)–target prior databases.
The panels in the left column show the boxplot of velocity consistency, the in-cluster coherence, and the cross-boundary direction correctness, respectively. The panels in the right column show the barplot of their mean values. For these boxplots in the left column, the whiskers extend to 1.5× the interquartile range of the distribution from the hinge. 3696 samples are included in the boxplot for each method.
Figure 1—figure supplement 5
The quantitative comparison of TSvelo modeling results on the gastrulation erythroid dataset using different transcription factor (TF)–target prior databases.
The panels in the left column show the boxplot of velocity consistency, the in-cluster coherence, and the cross-boundary direction correctness, respectively. The panels in the right column show the barplot of their mean values. For these boxplots in the left column, the whiskers extend to 1.5× the interquartile range of the distribution from the hinge. 9815 samples are included in the boxplot for each method.
Figure 1—figure supplement 6
Effect of the number of expectation–maximization (EM) iterations on TSvelo performance in the simulation dataset.
(a) Quantitative evaluation across different iteration settings. (b) 2D visualization of inferred dynamics. The results show that TSvelo converges rapidly, and performance remains stable after three iterations due to the built-in early stopping mechanism in the EM framework. For these boxplots in the top row, the whiskers extend to 1.5× the interquartile range of the distribution from the hinge. The sample size is 1800.
Figure 1—figure supplement 7
The predicted velocity streams of TSvelo, Dynamo, and cellDancer on the null simulation data.
Figure 2 with 3 supplements
Results on pancreas dataset.
(a) The pseudotime learned with TSvelo. (b) The stream plot for visualizing the RNA velocity inferred by TSvelo. (c) A quantitative comparison of cell-state separability between the 3D phase portrait used in TSvelo and the traditional 2D phase portrait. The down, central, and up hinges correspond to the first quartile, median value, and third quartile, respectively. The whiskers extend to 1.5× the interquartile range of the distribution from the hinge. 3696 samples are included in the boxplot for each method. (d) The dynamics fitting on Maml3, Anxa4, and Gstz1. For each gene, four plots are displayed in a 2 × 2 layout: the u–t, s–t, u–s, and alpha–u plots. (e) The unspliced–spliced phase portrait fitting on Maml3, Anxa4, and Gstz1 obtained by four baseline RNA velocity approaches, including scVelo, Dynamo, UniTVelo, and cellDancer. (f) The dynamics fitting in transcription–unspliced–spliced 3D phase portrait for Maml3, Anxa4, and Gstz1.
-
Figure 2—source data 1
Source data of the comparison between the 2D and 3D phase portrait fitting.
- https://cdn.elifesciences.org/articles/108950/elife-108950-fig2-data1-v1.csv
Figure 2—figure supplement 1
The velocity consistency, in-cluster coherence, and cross-boundary direction correctness on pancreas dataset.
The whiskers extend to 1.5× the interquartile range of the distribution from the hinge. 3696 samples are included in the boxplot for each method. In each plot, methods are ranked in descending order of their mean values. Numbers at the bottom indicate the sample size for each metric. Significance is determined using a one-sided Mann–Whitney U test. *****, ****, ***, **, and * represent p < 0.00001, 0.00001 ≤ p < 0.0001, 0.0001 ≤ p < 0.001, 0.001 ≤ p < 0.01, and 0.01 ≤ p < 0.05, respectively.
Figure 2—figure supplement 2
Dynamics fitting of velocity genes on the pancreas dataset.
For each gene, four plots are displayed in a 2 × 2 layout: the u–t, s–t, u–s, and alpha–u plots.
Figure 2—figure supplement 3
Dynamics fitting for additional transcription factors (TFs) on the pancreas dataset.
TSvelo directly models the dynamic between transcription and spliced abundance on these genes. For each gene, four plots are displayed in a 2 × 2 layout: the u–t, s–t, u–s, and alpha–u plots.
Figure 3 with 2 supplements
Results on gastrulation erythroid dataset.
(a) The Gene Ontology (GO) terms which are mostly enriched in the selected velocity genes of TSvelo. (b) The pseudotime learned with TSvelo. (c) The stream plot for visualizing the RNA velocity inferred by TSvelo. The quantitative comparison between TSvelo and multiple baseline approaches in terms of velocity consistency (d), in-cluster coherence (e), and cross-boundary direction correctness (f). The down, central, and up hinges correspond to the first quartile, median value, and third quartile, respectively. The whiskers extend to 1.5× the interquartile range of the distribution from the hinge. 9815 samples are included in the boxplot for each method. (g) The dynamics fitting on Hsp90ab1 obtained by TSvelo. Four plots are displayed in a 2 × 2 layout: the u–t, s–t, u–s, and alpha–u plots, where u, a, and alpha mean the abundance of unspliced mRNA, the abundance of spliced mRNA, and the learned transcriptional representation, respectively. (h) The phase portrait fitting on Hsp90ab1 obtained by baseline approaches. (i) The dynamics fitting on Rps26 obtained by TSvelo. (j) The transcription factors (TFs) with the highest ranked weights as identified by TSvelo. (k) Weights for Klf1’s targets, with the highest absolute weight. (l) Temporal dynamics of Klf1 and its target genes with the highest weights along pseudotime, which includes Hba-x, Alas2, and Gypa.
-
Figure 3—source data 1
The comparison of velocity consistency across different methods.
- https://cdn.elifesciences.org/articles/108950/elife-108950-fig3-data1-v1.zip
-
Figure 3—source data 2
The comparison of in-cluster coherence across different methods.
- https://cdn.elifesciences.org/articles/108950/elife-108950-fig3-data2-v1.zip
-
Figure 3—source data 3
The comparison of cross-boundary direction correctness across different methods.
- https://cdn.elifesciences.org/articles/108950/elife-108950-fig3-data3-v1.zip
Figure 3—figure supplement 1
The evaluation of the separability of cell-type labels in both the 2D (unspliced–spliced) phase portrait and the 3D (α–unspliced–spliced) phase portrait for the gastrulation erythroid dataset.
Figure 3—figure supplement 2
Dynamics fitting of velocity genes on the gastrulation erythroid dataset.
For each gene, four plots are displayed in a 2 × 2 layout: the u–t, s–t, u–s, and alpha–u plots.
Figure 4 with 3 supplements
Results on mouse brain dataset.
(a) The velocity stream inferred by Multivelo. (b) The stream plot for visualization of the RNA velocity inferred by TSvelo. (c) The pseudotime learned with TSvelo. The dynamics fitting on gene Meis2 (d), Basp1 (e), and Msi2 (f) by both Multivelo and TSvelo. In each panel, the leftmost plot shows the phase portrait fitting of Multivelo, and the next two columns show TSvelo’s results. The rightmost plot shows the learned transcriptional rate (in green), unspliced abundance (in blue), and spliced abundance (in red) along the pseudotime. Since the transcriptional rate is calculated for each individual cell, we apply a Generalized Additive Model (GAM) to transcriptional representation across cells along the pseudotime and present GAM-fitting results to better visualize its trends in these plots.
Figure 4—figure supplement 1
Dynamics fitting for additional genes on the mouse brain dataset.
For each gene, four plots are displayed in a 2 × 2 layout: the u–t, s–t, u–s, and alpha–u plots.
Figure 4—figure supplement 2
The distribution of gene-wise Spearman correlation between the estimated transcription rate in TSvelo and the chromatin accessibility rate in Multivelo.
Figure 4—figure supplement 3
Results on pons dataset.
(a) The Gene Ontology (GO) terms which are mostly enriched in the selected velocity genes of TSvelo. (b) The Leiden clustering. (c) The U-to-S delay under different initial cluster choices. (d) The pseudotime learned with TSvelo. (e) The stream plot for visualization of the RNA velocity inferred by TSvelo. (f) Dynamics fitting for multiple genes on the pons dataset. For each gene, four plots are displayed in a 2 × 2 layout: the u–t, s–t, u–s, and alpha–u plots.
Figure 5 with 2 supplements
Results on the multi-lineage dentate gyrus dataset.
(a) The Gene Ontology (GO) terms enriched in the selected velocity genes of TSvelo. (b) The pseudotime learned with TSvelo. (c) The stream plot for visualizing the RNA velocity inferred by TSvelo. Three lineages are detected, which are Granule lineage, CA lineage, and glial lineage. (d) The velocity stream inferred by scVelo. (e) The velocity stream inferred by cellDancer. The dynamics modeling of TSvelo on three axonogenesis-related genes, Ank3 (f), Map1b (g), and Slc1a2 (h). In the leftmost plot of each panel, the lines represent the predicted spliced abundance across all lineages, with the color indicating the cell types most associated with each pseudotime point along the corresponding lineage. Additionally, expression data for each lineage are shown as translucent points. The remaining plots in each panel display the dynamics of the learned transcriptional rate (in green), unspliced abundance (in blue), and spliced abundance (in red) along pseudotime for each lineage. The transcriptional representation in these plots is also processed using GAM fitting.
Figure 5—figure supplement 1
Lineages segmentation on the dentate gyrus dataset.
(a) The Leiden clustering. (b) The cluster–cluster graph constructed with PAGA. (c) The diffusion pseudotime when choosing initial cluster 3 as the initial one. (d) The three branches detected when choosing initial cluster 3 as the initial one. (e) U-to-S delay under different initial cluster choices. (f) The velocity stream and pseudotime inferred by TSvelo on each lineage. (g) The velocity stream obtained with the combined velocities from each lineage. (h) The combined TSvelo pseudotime from each lineage.
Figure 5—figure supplement 2
Dynamics fitting for genes Ank3, Map1b, and Slc1a2 on each lineage of the dentate gyrus dataset.
For each gene, four plots are displayed in a 2 × 2 layout: the u–t, s–t, u–s, and alpha–u plots.
Figure 6
Results on the LARRY dataset.
(a) The Leiden clustering on LARRY. (b) The initial Leiden cluster detection in preprocessing. (c) The pseudotime learned with TSvelo. (d) The stream plot for visualizing the RNA velocity inferred by TSvelo. (e) The Gene Ontology (GO) terms enriched in the selected velocity genes of TSvelo. (f) The dynamics modeling of TSvelo on four genes related to neutrophil development, which are Pygl, Ms4a3, Clec12a, and Lta4h. The lines represent the predicted spliced abundance across all lineages, with the color indicating the cell types most strongly associated with each pseudotime point along the corresponding lineage. (g) The dynamics modeling of Pygl, Ms4a3, Clec12a, and Lta4h on the neutrophil lineage.
Tables
Key resources table
| Reagent type (species) or resource | Designation | Source or reference | Identifiers | Additional information |
|---|---|---|---|---|
| Software, algorithm | TSvelo | https://github.com/lijc0804/TSvelo, copy archived at Li, 2026 |
Appendix 1—table 1
The comparison of RNA velocity approaches.
| Methods | Biological scope | Key assumptions | Inference framework | Modeling granularity |
|---|---|---|---|---|
| velocyto | Transcription and splicing | Steady state. | Least squares regression | Per-gene |
| scVelo | Transcription and splicing | Two-step transcription rates. Constant splicing and degradation rates. | Expectation–maximization (EM) | Per-gene |
| VeloAE | Transcription and splicing | Steady state. | Auto-encoder (AE) | Joint (all genes simultaneously) |
| UniTVelo | Transcription and splicing | Spliced RNA abundance is a twice-differentiable function of time. | Radial basis function (RBF) | Per-gene |
| Dynamo | Transcription and splicing or metabolic | Two-step transcription rates. Constant splicing and degradation rates. | Ordinary differential equations (ODEs) and sparseVFC | Per-gene |
| cellDancer | Transcription and splicing | Cell-specific time-dependent transcription, splicing, and degradation rates. | Convolutional neural network (CNN) | Per-gene |
| MultiVelo | Transcription and splicing | Two-step chromatin accessibility states. Constant splicing and degradation rates. | EM | Per-gene |
| BayVel | Transcription and splicing | A product of Poisson distributions for the initial state. | Markov chain Monte Carlo | Per-gene |
| TFvelo | Regulation | Cell-specific time-dependent transcription rate. Constant degradation rates. | Generalized expectation–maximization | Per-gene |
| TSvelo | Regulation, transcription, and splicing | Cell-specific time-dependent transcription rates. Constant splicing and degradation rates. | Neural ODE | Joint (all genes simultaneously) |
Appendix 1—table 2
Comparison of memory and runtime cost of different methods.
| Dataset scale | 600 cells | 1200 cells | 1800 cells | |
|---|---|---|---|---|
| scVelo | GPU memory | - | - | - |
| Modeling time | 10 s | 12 s | 14 s | |
| Dynamo | GPU memory | - | - | - |
| Modeling time | 24 s | 20 s | 24 s | |
| cellDancer | GPU memory | - | - | - |
| Modeling time | 82 s | 118 s | 245 s | |
| UniTVelo | GPU memory | 930 MB | 994 MB | 1122 MB |
| Modeling time | 345 s | 637 s | 1016 s | |
| TSvelo | GPU memory | 1262 MB | 1264 MB | 1264 MB |
| Modeling time | 734 s | 1544 s | 3027 s |
Appendix 1—table 3
Runtime of TSvelo under different numbers of expectation–maximization (EM) iterations on the simulation dataset.
The results demonstrate that TSvelo achieves stable performance even with a small number of iterations. Due to the early stopping mechanism of the EM framework, since TSvelo gets the best results at iteration 3, it actually performs no more than six iterations.
| Number of max iterations set in hyperparameters | 1 | 2 | 3 | 4 | 5 |
| Running time (s) | 1007 | 1168 | 2145 | 2208 | 2989 |
| Number of max iterations set in hyperparameters | 6 | 7 | 8 | 9 | 10 |
| Running time (s) | 3120 | 3256 | 3126 | 3121 | 3027 |
Additional files
Download links
A two-part list of links to download the article, or parts of the article, in various formats.
Downloads (link to download the article as PDF)
Open citations (links to open the citations from this article in various online reference manager services)
Cite this article (links to download the citations from this article in formats compatible with various reference manager tools)
Comprehensive RNA velocity by modeling the cascade of gene regulation, transcription, and splicing from single-cell RNA sequencing data with TSvelo
eLife 14:RP108950.
https://doi.org/10.7554/eLife.108950.4
