Kernel tricks

This example shows some niche, but useful functionalities of cellrank.tl.kernels.Kernel.

CellRank is split into cellrank.tl.kernels and :mod:`cellrank.tl.estimators:. Kernels compute transition matrices based on some inputs, like RNA velocity [Bergen20], [Manno18], while estimators perform inference based on a given kernel, e.g. they compute initial and terminal cells and fate probabilities.

Here, will will dive a bit deeper into the how these kernel objects work.

import cellrank as cr

adata = cr.datasets.pancreas_preprocessed("../example.h5ad")
adata

Out:

AnnData object with n_obs × n_vars = 2531 × 2000
    obs: 'day', 'proliferation', 'G2M_score', 'S_score', 'phase', 'clusters_coarse', 'clusters', 'clusters_fine', 'louvain_Alpha', 'louvain_Beta', 'initial_size_unspliced', 'initial_size_spliced', 'initial_size', 'n_counts', 'velocity_self_transition', 'dpt_pseudotime'
    var: 'highly_variable_genes', 'gene_count_corr', 'means', 'dispersions', 'dispersions_norm', 'fit_r2', 'fit_alpha', 'fit_beta', 'fit_gamma', 'fit_t_', 'fit_scaling', 'fit_std_u', 'fit_std_s', 'fit_likelihood', 'fit_u0', 'fit_s0', 'fit_pval_steady', 'fit_steady_u', 'fit_steady_s', 'fit_variance', 'fit_alignment_scaling', 'velocity_genes'
    uns: 'clusters_colors', 'clusters_fine_colors', 'diffmap_evals', 'iroot', 'louvain_Alpha_colors', 'louvain_Beta_colors', 'neighbors', 'pca', 'recover_dynamics', 'velocity_graph', 'velocity_graph_neg', 'velocity_params'
    obsm: 'X_diffmap', 'X_pca', 'X_umap', 'velocity_umap'
    varm: 'PCs', 'loss'
    layers: 'Ms', 'Mu', 'fit_t', 'fit_tau', 'fit_tau_', 'spliced', 'unspliced', 'velocity', 'velocity_u'
    obsp: 'connectivities', 'distances'

First, we create some kernels which will be used to compute the cell-to-cell transition matrix:

ck = cr.tl.kernels.ConnectivityKernel(adata)
vk = cr.tl.kernels.VelocityKernel(adata)
pk = cr.tl.kernels.PalantirKernel(adata)

Kernels can be easily combined with each other. All of the constants within parentheses (if they are needed) will be normalized to 1 automatically. This effect will be shown after the transition matrix has been computed.

10 * vk + 15 * ck

Out:

((10 * <Velo>) + (15 * <Conn>))

We can also build more complex kernel expression. Note that multiple usage of the same kernel instance in the expression caches it’s transition matrix, so it’s only computed once.

k = ((vk + vk + 42 * vk) + ck * 2) + ck * (3 * vk + 4 * pk)
k

Out:

((1 * ((1 * <Velo>) + (1 * <Velo>) + (42 * <Velo>) + (2 * <Conn>))) + ((1 * <Conn>) * (1 * ((3 * <Velo>) + (4 * <Pala>)))))

For complex expression like the one above, one can get the unique base kernels as follows.

k.kernels

Out:

[<Velo>, <Pala>, <Conn>]

Kernels can also be inverted (i.e. the backward attribute is set to the opposite value to the current one). Note that this operation does NOT create a copy and modifies the kernel inplace, most importantly, the computed transition matrix is removed. This makes it more easier to recompute the transition matrix of kernels, since the handles point to the original objects.

Note that in the 2nd print() statement, we access the private attribute - that’s because accessing cellrank.tl.kernels.Kernel.transition_matrix computes the transition matrix with default values. This happens only with basic kernels and not the kernel expressions and only if they are not part of a larger expression.

ck.compute_transition_matrix()
print(ck.transition_matrix is not None)

inv_ck = ~ck
print(ck._transition_matrix is not None)
print(inv_ck is ck)

Out:

True
False
True

We can also easily copy the kernel objects.

vk.compute_transition_matrix(
    mode="deterministic", softmax_scale=4, show_progress_bar=False
)
print(vk.transition_matrix is not None)

inv_vk = ~(vk.copy())
print(vk._transition_matrix is not None)
print(vk is inv_vk)

Out:

True
True
False

After the transition matrix is computed, we can inspect the parameters used in its computation.

ck.compute_transition_matrix()
ck.params

Out:

{'dnorm': True}

The transition matrix can be written to anndata.AnnData object. It is saved in the .obsp attribute. This also writes the parameters to the .uns attribute.

ck.write_to_adata(key="transition_matrix")
adata.obsp["transition_matrix"]

Out:

<2531x2531 sparse matrix of type '<class 'numpy.float64'>'
    with 102924 stored elements in Compressed Sparse Row format>

Precomputed kernels are useful if you have a transition matrix that was computed outside of CellRank that you would like to analyze using one of our estimators. It’s essentially an interface to CellRank’s estimators.

Below we supply a transition matrix saved in adata.obsp["transition_matrix"], but anndata.AnnData object is not required and we can just supply either numpy or scipy.sparse array. In that case, a minimal empty anndata.AnnData object is created, as shown below.

pk = cr.tl.kernels.PrecomputedKernel(adata.obsp["transition_matrix"])
pk.adata

Out:

AnnData object with n_obs × n_vars = 2531 × 1

Total running time of the script: ( 0 minutes 12.676 seconds)

Estimated memory usage: 219 MB

Gallery generated by Sphinx-Gallery