#!/usr/bin/env python
# coding: utf-8

# # DPT for hematopoiesis in mouse [(Moignard *et al.*, 2015)](http://dx.doi.org/10.1038/nbt.3154)

# The following analysis has been published by [Haghverdi *et al.*, Nat. Meth. (2016)](http://dx.doi.org/10.1038/nmeth.3971) together with a Matlab implementation of Diffusion Pseudotime. The data is from [Moignard *et al.*, Nat. Biotechn. (2015)](http://dx.doi.org/10.1038/nbt.3154).

# In[1]:


import numpy as np
import pandas as pd
import scanpy.api as sc

sc.settings.verbosity = 3  # verbosity: errors (0), warnings (1), info (2), hints (3)
sc.logging.print_versions()
results_file = './write/moignard15.h5ad'


# In[2]:


sc.settings.set_figure_params(dpi=80)  # low dpi (dots per inch) yields small inline figures


# This is single-cell qPCR data. It's already on a logrithmic scale and pretty low-dimensional. We do not have to perform extensive preprocessing.

# 1. Filter out a few genes.
# 2. Choose a 'root cell'.                                                                                                                                                           
# 3. Define groupnames by inspecting cellnames.                                                                                                                                    

# In[3]:


adata = sc.datasets.moignard15()


# Compute the neighborhood relations of single cells.

# In[4]:


sc.pp.neighbors(adata, n_neighbors=5, method='gauss', knn=False)


# Compute branchings and diffusion pseudotime using DPT.

# In[5]:


sc.tl.diffmap(adata)


# In[8]:


sc.tl.dpt(adata, n_branchings=1, n_dcs=10)


# In[9]:


sc.pl.diffmap(adata, color=['dpt_pseudotime', 'dpt_groups', 'exp_groups'])


# In[10]:


sc.pl.dpt_groups_pseudotime(adata)


# In[11]:


sc.pl.dpt_timeseries(adata)


# Let us annotate the cell groups as follows.

# In[13]:


adata.rename_categories('dpt_groups', ['undecided/endothelial', 'erythrocytes', 'trunk', 'endothelial'])


# In[14]:


sc.pl.diffmap(adata, color='dpt_groups')


# Save the results.

# In[15]:


adata.write(results_file)


# Or as csv.

# In[11]:


# adata.obs.to_csv('./write/annotation.csv')


# In[12]:


# To write the full object to csvs
# adata.write_csvs('./write/moignard15.csv')