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

# # Example - FRET histogram fitting
# 
# *This notebook is part of smFRET burst analysis software [FRETBursts](http://opensmfs.github.io/FRETBursts/).*

# > In this notebook shows how to fit a FRET histogram.
# > For a complete tutorial on burst analysis see 
# > [FRETBursts - us-ALEX smFRET burst analysis](FRETBursts - us-ALEX smFRET burst analysis.ipynb).

# In[1]:


from fretbursts import *


# In[2]:


sns = init_notebook(apionly=True)
import lmfit
print('lmfit version:', lmfit.__version__)


# In[3]:


# Tweak here matplotlib style
import matplotlib as mpl
mpl.rcParams['font.sans-serif'].insert(0, 'Arial')
mpl.rcParams['font.size'] = 12
get_ipython().run_line_magic('config', "InlineBackend.figure_format = 'retina'")


# # Get and process data

# In[4]:


url = 'http://files.figshare.com/2182601/0023uLRpitc_NTP_20dT_0.5GndCl.hdf5'
download_file(url, save_dir='./data')
full_fname = "./data/0023uLRpitc_NTP_20dT_0.5GndCl.hdf5"

d = loader.photon_hdf5(full_fname)
loader.alex_apply_period(d)
d.calc_bg(bg.exp_fit, time_s=1000, tail_min_us=(800, 4000, 1500, 1000, 3000))
d.burst_search(L=10, m=10, F=6)
ds = d.select_bursts(select_bursts.size, add_naa=True, th1=30)


# # Fitting the FRET histogram
# 
# 

# We start defining the model. Here we choose a 3-Gaussian model:

# In[5]:


model = mfit.factory_three_gaussians()
model.print_param_hints()


# The previsou cell prints all the model parameters. 
# Each parameters has an initial value and bounds (min, max).
# The column `vary` tells if a parameter is varied during the fit
# (if False the parameter is fixed).
# Parameters with an expression (`Expr` column) are not free but
# the are computed as a function of other parameters.
# 
# We can modify the paramenters constrains as follows:

# In[6]:


model.set_param_hint('p1_center', value=0.1, min=-0.1, max=0.3)
model.set_param_hint('p2_center', value=0.4, min=0.3, max=0.7)
model.set_param_hint('p2_sigma', value=0.04, min=0.02, max=0.18)
model.set_param_hint('p3_center', value=0.85, min=0.7, max=1.1)


# Then, we fit and plot the model:

# In[7]:


E_fitter = bext.bursts_fitter(ds, 'E', binwidth=0.03)
E_fitter.fit_histogram(model=model, pdf=False, method='nelder')
E_fitter.fit_histogram(model=model, pdf=False, method='leastsq')
dplot(ds, hist_fret, show_model=True, pdf=False);


# In[8]:


# dplot(ds, hist_fret, show_model=True, pdf=False, figsize=(6, 4.5));
# plt.xlim(-0.1, 1.1)
# plt.savefig('fret_hist_fit.png', bbox_inches='tight', dpi=200, transparent=False)


# The results are in `E_fitter`:

# In[9]:


res = E_fitter.fit_res[0]
res.params.pretty_print()


# To get a dictionary of values:

# In[10]:


res.values


# This is the startndard lmfit's fit report:

# In[11]:


print(res.fit_report(min_correl=0.5))


# The previous cell reports error ranges computed from the covariance matrix.
# [More accurare confidence intervals](https://lmfit.github.io/lmfit-py/confidence.html) 
# can be obtained with:

# In[12]:


ci = res.conf_interval()
lmfit.report_ci(ci)


# # Tidy fit results
# 
# It is convenient to put the fit results in a DataFrame for further analysis.
# A dataframe of fitted parameters is already in `E_fitter`:

# In[13]:


E_fitter.params


# With [`pybroom`](http://pybroom.readthedocs.io/) we can get a "tidy" DataFrame
# with more complete fit results:

# In[14]:


import pybroom as br


# In[15]:


df = br.tidy(res)
df


# Now, for example, we can easily select parameters by name:

# In[16]:


df.loc[df.name.str.contains('center')]


# In[17]:


df.loc[df.name.str.contains('sigma')]