Fixed point analysis - Chase and pulse experiment#
Dynamo’s fixed-point-based analysis of scEU-seq organoid data using only the chase and pulse experiment.
Library imports#
import sys
import pandas as pd
import dynamo as dyn
import scanpy as sc
import scvelo as scv
sys.path.extend(["../../../", "."])
from paths import DATA_DIR # isort: skip # noqa: E402
|-----> setting visualization default mode in dynamo. Your customized matplotlib settings might be overritten.
/vol/storage/miniconda3/envs/dynamo-py39/lib/python3.9/site-packages/nxviz/__init__.py:18: UserWarning:
nxviz has a new API! Version 0.7.3 onwards, the old class-based API is being
deprecated in favour of a new API focused on advancing a grammar of network
graphics. If your plotting code depends on the old API, please consider
pinning nxviz at version 0.7.3, as the new API will break your old code.
To check out the new API, please head over to the docs at
https://ericmjl.github.io/nxviz/ to learn more. We hope you enjoy using it!
(This deprecation message will go away in version 1.0.)
warnings.warn(
Global seed set to 0
General settings#
sc.settings.verbosity = 3
scv.settings.verbosity = 3
Data loading#
adata = sc.read(DATA_DIR / "sceu_organoid" / "processed" / "adata_dynamo-chase_and_pulse-2000features.h5ad")
adata
AnnData object with n_obs × n_vars = 3452 × 2000
obs: 'experiment', 'time', 'cell_type', 'som_cluster_id', 'cell_type_merged', 'initial_size', 'n_counts', 'pass_basic_filter', 'ntr'
var: 'ensum_id', 'gene_count_corr', 'means', 'dispersions', 'dispersions_norm', 'highly_variable', 'pass_basic_filter', 'use_for_pca', 'ntr', 'alpha', 'a', 'b', 'alpha_a', 'alpha_i', 'beta', 'p_half_life', 'gamma', 'half_life', 'cost', 'logLL', 'use_for_dynamics'
uns: 'PCs', 'cell_type_colors', 'dynamics', 'neighbors', 'pca', 'pca_fit', 'pca_mean', 'pp', 'umap'
obsm: 'X_pca', 'X_umap', 'X_umap_paper'
varm: 'PCs'
layers: 'M_n', 'M_nn', 'M_t', 'M_tn', 'M_tt', 'X_new', 'X_total', 'new', 'total', 'velocity_N', 'velocity_T'
obsp: 'connectivities', 'distances', 'moments_con'
Preprocessing#
dyn.tl.reduceDimension(adata, layer="X_new", enforce=True)
|-----> retrive data for non-linear dimension reduction...
|-----> perform umap...
|-----> [dimension_reduction projection] in progress: 100.0000%
|-----> [dimension_reduction projection] finished [24.1703s]
Fixed point analysis#
Velocity with total RNA#
dyn.tl.cell_velocities(adata, ekey="M_t", vkey="velocity_T", enforce=True)
dyn.vf.VectorField(adata, basis="umap")
dyn.vf.topography(adata)
|-----> 0 genes are removed because of nan velocity values.
|-----> [calculating transition matrix via pearson kernel with sqrt transform.] in progress: 100.0000%
|-----> [calculating transition matrix via pearson kernel with sqrt transform.] finished [1.9442s]
|-----> [projecting velocity vector to low dimensional embedding] in progress: 100.0000%
|-----> [projecting velocity vector to low dimensional embedding] finished [0.8892s]
|-----> VectorField reconstruction begins...
|-----> Retrieve X and V based on basis: UMAP.
Vector field will be learned in the UMAP space.
|-----> Generating high dimensional grids and convert into a row matrix.
|-----> Learning vector field with method: sparsevfc.
|-----> [SparseVFC] begins...
|-----> Sampling control points based on data velocity magnitude...
|-----> [SparseVFC] in progress: 100.0000%
|-----> [SparseVFC] finished [0.5633s]
|-----> <insert> velocity_umap_SparseVFC to obsm in AnnData Object.
|-----> <insert> X_umap_SparseVFC to obsm in AnnData Object.
|-----> <insert> VecFld_umap to uns in AnnData Object.
|-----> <insert> control_point_umap to obs in AnnData Object.
|-----> <insert> inlier_prob_umap to obs in AnnData Object.
|-----> <insert> obs_vf_angle_umap to obs in AnnData Object.
|-----> [VectorField] in progress: 100.0000%
|-----> [VectorField] finished [0.7333s]
AnnData object with n_obs × n_vars = 3452 × 2000
obs: 'experiment', 'time', 'cell_type', 'som_cluster_id', 'cell_type_merged', 'initial_size', 'n_counts', 'pass_basic_filter', 'ntr', 'control_point_umap', 'inlier_prob_umap', 'obs_vf_angle_umap'
var: 'ensum_id', 'gene_count_corr', 'means', 'dispersions', 'dispersions_norm', 'highly_variable', 'pass_basic_filter', 'use_for_pca', 'ntr', 'alpha', 'a', 'b', 'alpha_a', 'alpha_i', 'beta', 'p_half_life', 'gamma', 'half_life', 'cost', 'logLL', 'use_for_dynamics', 'use_for_transition'
uns: 'PCs', 'cell_type_colors', 'dynamics', 'neighbors', 'pca', 'pca_fit', 'pca_mean', 'pp', 'umap', 'explained_variance_ratio_', 'pca_valid_ind', 'X_new_neighbors', 'umap_fit', 'grid_velocity_umap', 'VecFld_umap'
obsm: 'X_pca', 'X_umap', 'X_umap_paper', 'X_new_pca', 'X_new_umap', 'velocity_umap', 'velocity_umap_SparseVFC', 'X_umap_SparseVFC'
varm: 'PCs'
layers: 'M_n', 'M_nn', 'M_t', 'M_tn', 'M_tt', 'X_new', 'X_total', 'new', 'total', 'velocity_N', 'velocity_T'
obsp: 'connectivities', 'distances', 'moments_con', 'pearson_transition_matrix'
df = pd.DataFrame(adata.uns["VecFld_umap"]["Xss"], columns=["umap_1", "umap_2"])
df["fixed_point_type"] = adata.uns["VecFld_umap"]["ftype"]
df["fixed_point_type"].replace({-1: "stable", 0: "saddle", 1: "unstable"}, inplace=True)
neighbor_idx = dyn.tools.utils.nearest_neighbors(df[["umap_1", "umap_2"]], adata.obsm["X_umap"])
df["cell_type"] = [
adata.obs.loc[adata.obs_names[neighbors], "cell_type"].mode().values[0] for neighbors in neighbor_idx
]
print(f"Stable FPs: {df.loc[df['fixed_point_type'] == 'stable', 'cell_type'].unique()}")
Stable FPs: ['Stem cells']
Velocity with new RNA#
dyn.tl.cell_velocities(adata, ekey="M_n", vkey="velocity_N", enforce=True)
dyn.vf.VectorField(adata, basis="umap")
dyn.vf.topography(adata)
|-----> retrive data for non-linear dimension reduction...
|-----? adata already have basis umap. dimension reduction umap will be skipped!
set enforce=True to re-performing dimension reduction.
|-----> [dimension_reduction projection] in progress: 100.0000%
|-----> [dimension_reduction projection] finished [0.0019s]
|-----> 0 genes are removed because of nan velocity values.
|-----> [calculating transition matrix via pearson kernel with sqrt transform.] in progress: 100.0000%
|-----> [calculating transition matrix via pearson kernel with sqrt transform.] finished [2.3528s]
|-----> [projecting velocity vector to low dimensional embedding] in progress: 100.0000%
|-----> [projecting velocity vector to low dimensional embedding] finished [0.8721s]
|-----> VectorField reconstruction begins...
|-----> Retrieve X and V based on basis: UMAP.
Vector field will be learned in the UMAP space.
|-----> Generating high dimensional grids and convert into a row matrix.
|-----> Learning vector field with method: sparsevfc.
|-----> [SparseVFC] begins...
|-----> Sampling control points based on data velocity magnitude...
|-----> [SparseVFC] in progress: 100.0000%
|-----> [SparseVFC] finished [0.5459s]
|-----------> current cosine correlation between input velocities and learned velocities is less than 0.6. Make a 1-th vector field reconstruction trial.
|-----> [SparseVFC] begins...
|-----> Sampling control points based on data velocity magnitude...
|-----> [SparseVFC] in progress: 100.0000%
|-----> [SparseVFC] finished [0.7590s]
|-----------> current cosine correlation between input velocities and learned velocities is less than 0.6. Make a 2-th vector field reconstruction trial.
|-----> [SparseVFC] begins...
|-----> Sampling control points based on data velocity magnitude...
|-----> [SparseVFC] in progress: 100.0000%
|-----> [SparseVFC] finished [0.6202s]
|-----------> current cosine correlation between input velocities and learned velocities is less than 0.6. Make a 3-th vector field reconstruction trial.
|-----> [SparseVFC] begins...
|-----> Sampling control points based on data velocity magnitude...
|-----> [SparseVFC] in progress: 100.0000%
|-----> [SparseVFC] finished [0.6246s]
|-----------> current cosine correlation between input velocities and learned velocities is less than 0.6. Make a 4-th vector field reconstruction trial.
|-----> [SparseVFC] begins...
|-----> Sampling control points based on data velocity magnitude...
|-----> [SparseVFC] in progress: 100.0000%
|-----> [SparseVFC] finished [0.6225s]
|-----------> current cosine correlation between input velocities and learned velocities is less than 0.6. Make a 5-th vector field reconstruction trial.
|-----? Cosine correlation between input velocities and learned velocities is less than 0.6 after 5 trials of vector field reconstruction.
|-----> <insert> velocity_umap_SparseVFC to obsm in AnnData Object.
|-----> <insert> X_umap_SparseVFC to obsm in AnnData Object.
|-----> <insert> VecFld_umap to uns in AnnData Object.
|-----> <insert> control_point_umap to obs in AnnData Object.
|-----> <insert> inlier_prob_umap to obs in AnnData Object.
|-----> <insert> obs_vf_angle_umap to obs in AnnData Object.
|-----> [VectorField] in progress: 100.0000%
|-----> [VectorField] finished [3.4083s]
AnnData object with n_obs × n_vars = 3452 × 2000
obs: 'experiment', 'time', 'cell_type', 'som_cluster_id', 'cell_type_merged', 'initial_size', 'n_counts', 'pass_basic_filter', 'ntr', 'control_point_umap', 'inlier_prob_umap', 'obs_vf_angle_umap'
var: 'ensum_id', 'gene_count_corr', 'means', 'dispersions', 'dispersions_norm', 'highly_variable', 'pass_basic_filter', 'use_for_pca', 'ntr', 'alpha', 'a', 'b', 'alpha_a', 'alpha_i', 'beta', 'p_half_life', 'gamma', 'half_life', 'cost', 'logLL', 'use_for_dynamics', 'use_for_transition'
uns: 'PCs', 'cell_type_colors', 'dynamics', 'neighbors', 'pca', 'pca_fit', 'pca_mean', 'pp', 'umap', 'explained_variance_ratio_', 'pca_valid_ind', 'X_new_neighbors', 'umap_fit', 'grid_velocity_umap', 'VecFld_umap'
obsm: 'X_pca', 'X_umap', 'X_umap_paper', 'X_new_pca', 'X_new_umap', 'velocity_umap', 'velocity_umap_SparseVFC', 'X_umap_SparseVFC'
varm: 'PCs'
layers: 'M_n', 'M_nn', 'M_t', 'M_tn', 'M_tt', 'X_new', 'X_total', 'new', 'total', 'velocity_N', 'velocity_T'
obsp: 'connectivities', 'distances', 'moments_con', 'pearson_transition_matrix'
df = pd.DataFrame(adata.uns["VecFld_umap"]["Xss"], columns=["umap_1", "umap_2"])
df["fixed_point_type"] = adata.uns["VecFld_umap"]["ftype"]
df["fixed_point_type"].replace({-1: "stable", 0: "saddle", 1: "unstable"}, inplace=True)
neighbor_idx = dyn.tools.utils.nearest_neighbors(df[["umap_1", "umap_2"]], adata.obsm["X_new_umap"])
df["cell_type"] = [
adata.obs.loc[adata.obs_names[neighbors], "cell_type"].mode().values[0] for neighbors in neighbor_idx
]
print(f"Stable FPs: {df.loc[df['fixed_point_type'] == 'stable', 'cell_type'].unique()}")
Stable FPs: ['Enterocytes']