PyBroom Example - Simple

This notebook is part of pybroom.

This notebook shows the simplest usage of pybroom when performing a curve fit of a single dataset. Possible applications are only hinted. For a more complex (and interesting!) example using multiple datasets see pybroom-example-multi-datasets.



In [ ]:

    
import numpy as np
from numpy import sqrt, pi, exp, linspace
from lmfit import Model
import matplotlib.pyplot as plt
%matplotlib inline
%config InlineBackend.figure_format='retina'  # for hi-dpi displays



In [ ]:

    
import lmfit
print('lmfit: %s' % lmfit.__version__)



In [ ]:

    
import pybroom as br

Create Noisy Data



In [ ]:

    
x = np.linspace(-10, 10, 101)



In [ ]:

    
peak1 = lmfit.models.GaussianModel(prefix='p1_')
peak2 = lmfit.models.GaussianModel(prefix='p2_')
model = peak1 + peak2



In [ ]:

    
params = model.make_params(p1_amplitude=1, p2_amplitude=1, 
                           p1_sigma=1, p2_sigma=1)



In [ ]:

    
y_data = model.eval(x=x, p1_center=-1, p2_center=2, p1_sigma=0.5, p2_sigma=1, p1_amplitude=1, p2_amplitude=2)
y_data.shape



In [ ]:

    
y_data += np.random.randn(*y_data.shape)/10



In [ ]:

    
plt.plot(x, y_data)

Model Fitting



In [ ]:

    
params = model.make_params(p1_center=0, p2_center=3, 
                           p1_sigma=0.5, p2_sigma=1, 
                           p1_amplitude=1, p2_amplitude=2)
result = model.fit(y_data, x=x, params=params)

Fit result from an lmfit Model can be inspected with with fit_report or params.pretty_print():



In [ ]:

    
print(result.fit_report())



In [ ]:

    
result.params.pretty_print()

These methods a re convenient but extracting the data from the lmfit object requires some work and the knowledge of lmfit object structure.

pybroom comes to help, extracting data from fit results and returning pandas DataFrame in tidy format that can be much more easily manipulated, filtered and plotted.

Glance

Glancing at the fit results (dropping some verbose columns):



In [ ]:

    
dg = br.glance(result)
dg.drop('model', 1).drop('message', 1)

The glance function returns a DataFrame with one row per fit-result object.

Application Idea

If you fit N models to the same dataset you can compare statistics such as reduced-$\chi^2$

Or, fitting several with several methods (and datasets) you can study the convergence properties using reduced-$\chi^2$, number of function evaluation and success rate.

Tidy

Tidy fit results for all the parameters:



In [ ]:

    
dt = br.tidy(result)
dt

The tidy function returns one row for each parameter.



In [ ]:

    
dt.loc[dt.name == 'p1_center']

Augment

Tidy dataframe with data function of the independent variable ('x'). Columns include the data being fitted, best fit, best fit components, residuals, etc.



In [ ]:

    
da = br.augment(result)
da.head()

The augment function returns one row for each data point.



In [ ]:

    
d = br.augment(result)



In [ ]:

    
fig, ax = plt.subplots(2, 1, figsize=(7, 8))
ax[1].plot('x', 'data', data=d, marker='o', ls='None')
ax[1].plot('x', "Model(gaussian, prefix='p1_')", data=d, lw=2, ls='--')
ax[1].plot('x', "Model(gaussian, prefix='p2_')", data=d, lw=2, ls='--')
ax[1].plot('x', 'best_fit', data=d, lw=2)
ax[0].plot('x', 'residual', data=d);

Application Idea

Fitting N datasets with the same model (or N models with the same dataset) you can automatically build a panel plot with seaborn using the dataset (or the model) as categorical variable. This example is illustrated in pybroom-example-multi-datasets.