Python implementation of the CIARA algorithm that integrates into scanpy's analysis with the AnnData format.
The package can be installed via pip:
python -m pip install ciara_python
Note: The package only works on UNIX / MacOS operating systems, not on Windows systems, due to the copy-on-write multiprocessing setup used.
The functions are designed to work with scanpy's AnnData objects. For an interactive tutorial check out the Human Gastrula IPython Notebook which is also part of this repository.
First you load your dataset as a scanpy object and after normal preprocessing calculate the knn-network:
import scanpy as sc
pbmc = sc.datasets.pbmc3k()
sc.pp.filter_cells(pbmc, min_genes=50)
sc.pp.filter_genes(pbmc, min_cells=10)
sc.pp.log1p(pbmc)
sc.pp.pca(pbmc)
sc.pp.neighbors(pbmc)
The CIARA package contains the two main functions get_background_full()
and ciara()
which should be imported via:
from ciara_python import get_background_full, ciara
Afterwards, the background genes get marked by running the get_background_full()
function on your scanpy dataset. This adds the boolean column 'CIARA_background' to your pbmc.var
AnnData slot, where relevant background genes are marked.
get_background_full(pbmc, threshold=1, n_cells_low=2, n_cells_high=5)
Finally, the ciara()
function is run on the dataset. This adds the column 'CIARA_p_value' to your pbmc.var
object, where the calculated p_values for each of the previously marked background genes are stored.
ciara(pbmc, n_cores=4, p_value=0.001, approximation=True, local_region=1)
We can then extract the top markers or rare cells (lowest CIARA_p_value) and plot them onto the UMAP:
from matplotlib.pyplot import rc_context
sc.tl.umap(pbmc)
top_markers = pbmc.var.nsmallest(4, ["CIARA_p_value"],)
with rc_context({'figure.figsize': (3, 3)}):
sc.pl.umap(pbmc, color=top_markers.index.tolist())
Link to R package: