A question about ‘g_fwd.compute_fate_probabilities(preconditioner='ilu', tol=1e-7, solver='direct')’

[liuyang218622]

Hi, when I run this code below:

import sys
import cellrank as cr
import scanpy as sc
import pandas as pd
sc.settings.set_figure_params(frameon=False, dpi=100)
import numpy as np
cr.settings.verbosity = 2
import warnings
warnings.simplefilter("ignore", category=UserWarning)

Load data

adata=sc.read_h5ad('data.h5ad')
sc.pl.embedding(adata, basis='velocity_umap', color=['celltype','leiden_velocity'],
legend_loc='on data')

Prepare Cellrank

sc.tl.pca(adata, n_comps=50, svd_solver = 'arpack')
sc.pp.neighbors(adata, n_neighbors=15, n_pcs=40)
sc.tl.diffmap(adata)
root = adata.obs.index.get_loc(adata[adata.obs['celltype']=='C13'].obs.index[np.argmin(adata[adata.obs['celltype']=='C13'].obsm['X_umap'][:, 1])])
adata.uns['iroot'] = root
sc.tl.dpt(adata)
adata.obs['root'] = np.nan
adata.obs['root'][root] = 1
sc.pl.embedding(adata, basis='velocity_umap', color=['dpt_pseudotime','celltype'],vmax=0.5)

Computing fate probability

from cellrank.kernels import PseudotimeKernel
pk = PseudotimeKernel(adata,time_key="dpt_pseudotime")
pk.compute_transition_matrix(threshold_scheme='soft',nu=0.01, b=20.0)
ck = cr.kernels.ConnectivityKernel(adata)
ck.compute_transition_matrix()
kernel = 0.3ck + 0.7pk
from cellrank.estimators import GPCCA
g_fwd = GPCCA(kernel)
print(g_fwd)
g_fwd.fit(n_states=3, cluster_key="celltype")

predict_terminal_status

g_fwd.plot_macrostates(which="all")
g_fwd.predict_terminal_states(method="top_n", n_states=5)
g_fwd.plot_macrostates(which="terminal")
CD4_T_cells = (adata[adata.obs['celltype'] == 'CD4 T cells']).obs_names[:30]
CD8_T_cells = (adata[adata.obs['celltype'] == 'CD8 T cells']).obs_names[:30]
Tgd_cells = (adata[adata.obs['celltype'] == 'Tgd cells']).obs_names[:30]
tLable = pd.DataFrame({
'CD4 T cells': CD4_T_cells,
'CD8 T cells': CD8_T_cells,
'Tgd cells': Tgd_cells
})
tLable = tLable.to_dict(orient='list')
g_fwd.set_terminal_states(states=tLable, cluster_key='celltype')
g_fwd.plot_macrostates(which="terminal")

Compute fate probabilities and driver genes

g_fwd.compute_fate_probabilities(preconditioner='ilu', tol=1e-7, solver='direct')

I got this error

Computing fate probabilities
100%
3/3 [19:59<00:00, 398.50s/]
WARNING: 3 solution(s) did not converge

ValueError Traceback (most recent call last)
Cell In[69], line 1
----> 1 g_fwd.compute_fate_probabilities(preconditioner='ilu', tol=1e-7, solver='direct')

File ~/anaconda3/envs/cellrank/lib/python3.10/site-packages/cellrank/estimators/mixins/_fate_probabilities.py:220, in FateProbsMixin.compute_fate_probabilities(self, keys, solver, use_petsc, n_jobs, backend, show_progress_bar, tol, preconditioner)
218 start = logg.info("Computing fate probabilities")
219 data = self._rec_trans_states(keys, ctx="fate_probs")
--> 220 fate_probs = self._compute_fate_probabilities(
221 data.q,
222 data.s,
223 trans_indices=data.trans_indices,
224 term_states=data.term_states,
225 solver=solver,
226 use_petsc=use_petsc,
227 n_jobs=n_jobs,
228 backend=backend,
229 tol=tol,
230 show_progress_bar=show_progress_bar,
231 preconditioner=preconditioner,
232 )
233 fate_probs = Lineage(
234 fate_probs,
235 names=list(data.term_states.cat.categories),
236 colors=data.term_states_colors,
237 )
239 params = self._create_params(remove=["use_petsc", "n_jobs", "backend", "show_progress_bar"])

File ~/anaconda3/envs/cellrank/lib/python3.10/site-packages/cellrank/estimators/mixins/_fate_probabilities.py:494, in FateProbsMixin._compute_fate_probabilities(self, q, s, trans_indices, term_states, solver, use_petsc, n_jobs, backend, tol, show_progress_bar, preconditioner)
492 mask = np.isclose(abs_classes.sum(1), 1.0, rtol=1e-3)
493 if not np.all(mask):
--> 494 raise ValueError(
495 f"{np.sum(~mask)} value(s) do not sum to 1 (rtol=1e-3). Try decreasing the tolerance as tol=..., "
496 f"specifying a preconditioner as preconditioner=... or "
497 f"use a direct solver as solver='direct' if the matrix is small."
498 )
500 return abs_classes

ValueError: 4791 value(s) do not sum to 1 (rtol=1e-3). Try decreasing the tolerance as tol=..., specifying a preconditioner as preconditioner=... or use a direct solver as solver='direct' if the matrix is small.

I try to set the tol=1e-7， tol=1e-20，tol=1e-25. I also used g_fwd = GPCCA(kernel) instead of g_fwd = GPCCA(pk). But I still get the error. Can you give me some suggestions?

[WeilerP]

@liuyang218622, can you check values larger than 1e-7? There was a counterintuitive way to decrease the threshold; I'm just not sure if it was here or if we updated it. Please also make sure you are running the latest version of CellRank.

[liuyang218622]

The cellrank version I used is v2.0.4

I used the tol=1e-2，tol=1e-3，tol=1e-5 but still get the same error

Computing fate probabilities
WARNING: Unable to import petsc4py. For installation, please refer to: https://petsc4py.readthedocs.io/en/stable/install.html.
Defaulting to 'gmres' solver.
100%
 3/3 [00:00<00:00, 15.53/s]

ValueError Traceback (most recent call last)
Cell In[29], line 1
----> 1 g_fwd.compute_fate_probabilities(preconditioner='ilu', tol=1e-2, solver='direct')

File ~/anaconda3/envs/Cellrank2/lib/python3.10/site-packages/cellrank/estimators/mixins/_fate_probabilities.py:220, in FateProbsMixin.compute_fate_probabilities(self, keys, solver, use_petsc, n_jobs, backend, show_progress_bar, tol, preconditioner)
218 start = logg.info("Computing fate probabilities")
219 data = self._rec_trans_states(keys, ctx="fate_probs")
--> 220 fate_probs = self._compute_fate_probabilities(
221 data.q,
222 data.s,
223 trans_indices=data.trans_indices,
224 term_states=data.term_states,
225 solver=solver,
226 use_petsc=use_petsc,
227 n_jobs=n_jobs,
228 backend=backend,
229 tol=tol,
230 show_progress_bar=show_progress_bar,
231 preconditioner=preconditioner,
232 )
233 fate_probs = Lineage(
234 fate_probs,
235 names=list(data.term_states.cat.categories),
236 colors=data.term_states_colors,
237 )
239 params = self._create_params(remove=["use_petsc", "n_jobs", "backend", "show_progress_bar"])

File ~/anaconda3/envs/Cellrank2/lib/python3.10/site-packages/cellrank/estimators/mixins/_fate_probabilities.py:494, in FateProbsMixin._compute_fate_probabilities(self, q, s, trans_indices, term_states, solver, use_petsc, n_jobs, backend, tol, show_progress_bar, preconditioner)
492 mask = np.isclose(abs_classes.sum(1), 1.0, rtol=1e-3)
493 if not np.all(mask):
--> 494 raise ValueError(
495 f"{np.sum(~mask)} value(s) do not sum to 1 (rtol=1e-3). Try decreasing the tolerance as tol=..., "
496 f"specifying a preconditioner as preconditioner=... or "
497 f"use a direct solver as solver='direct' if the matrix is small."
498 )
500 return abs_classes

ValueError: 6411 value(s) do not sum to 1 (rtol=1e-3). Try decreasing the tolerance as tol=..., specifying a preconditioner as preconditioner=... or use a direct solver as solver='direct' if the matrix is small.

[WeilerP]

@liuyang218622, could you please use CellRank with petsc and slepsc installed, and check again?

[liuyang218622]

I am installing them now, but I have another question. Is this condition associated with the fact that I set the terminal status manually?

Predict_terminal_status

g_fwd.plot_macrostates(which="all")
g_fwd.predict_terminal_states(method="top_n", n_states=5)
g_fwd.plot_macrostates(which="terminal")
CD4_T_cells = (adata[adata.obs['celltype'] == 'CD4 T cells']).obs_names[:30]
CD8_T_cells = (adata[adata.obs['celltype'] == 'CD8 T cells']).obs_names[:30]
Tgd_cells = (adata[adata.obs['celltype'] == 'Tgd cells']).obs_names[:30]
tLable = pd.DataFrame({
'CD4 T cells': CD4_T_cells,
'CD8 T cells': CD8_T_cells,
'Tgd cells': Tgd_cells
})
tLable = tLable.to_dict(orient='list')

g_fwd.set_terminal_states(states=tLable, cluster_key='celltype')

g_fwd.plot_macrostates(which="terminal")

[liuyang218622]

And I found If I don't set terminal_status manually. I can run g_fwd.compute_fate_probabilities(preconditioner='ilu', tol=1e-7, solver='direct') successfully