.. DO NOT EDIT. .. THIS FILE WAS AUTOMATICALLY GENERATED BY SPHINX-GALLERY. .. TO MAKE CHANGES, EDIT THE SOURCE PYTHON FILE: .. "auto_examples/fluorescence_decay/plot_fit23_benchmark_1.py" .. LINE NUMBERS ARE GIVEN BELOW. .. only:: html .. note:: :class: sphx-glr-download-link-note :ref:`Go to the end ` to download the full example code .. rst-class:: sphx-glr-example-title .. _sphx_glr_auto_examples_fluorescence_decay_plot_fit23_benchmark_1.py: =============================== MLE fit23: Benchmark =============================== This is a benchmark of ``fit2x.DecayFit.fit23`` that implements an maximum likelihood estimator to determine a fluorescence lifetime and a rotational correlation time for time and polarization resolved fluorescence decays with low photon counts. In this benchmark fluorescence decays are simulated for a fluorescence lifetime in the range ``(tau_min, tau_max)`` with varying number of photons. The simulated and the recovered fluorescence lifetimes are compared. .. GENERATED FROM PYTHON SOURCE LINES 13-137 .. image-sg:: /auto_examples/fluorescence_decay/images/sphx_glr_plot_fit23_benchmark_1_001.png :alt: Example decay, $\tau_{recov}$ for 20 photons, Photon count dependence :srcset: /auto_examples/fluorescence_decay/images/sphx_glr_plot_fit23_benchmark_1_001.png :class: sphx-glr-single-img .. code-block:: Python import tttrlib import numpy as np import pylab as plt np.random.seed(0) # setup some parameters tau_min = 0.1 tau_max = 5.0 n_photons_min = 20 n_photons_max = 120 n_photon_step = 2 n_samples = 200 n_channels = 32 irf_position_p = 2.0 irf_position_s = 18.0 irf_width = 0.25 period, g_factor, l1, l2, conv_stop = 32, 1.0, 0.1, 0.1, 31 dt = 0.5079365079365079 np.random.seed(0) irf_np = np.array([0, 0, 0, 260, 1582, 155, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 22, 1074, 830, 10, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0], dtype=np.float64) bg = np.zeros_like(irf_np) fit23 = tttrlib.Fit23( dt=dt, irf=irf_np, background=bg, period=period, g_factor=g_factor, l1=l1, l2=l2 ) tau, gamma, r0, rho = 2.0, 0.01, 0.38, 1.22 x0 = np.array([tau, gamma, r0, rho]) fixed = np.array([0, 1, 1, 0]) """ In this loop the fluorescence decays are simulated and the simulated decays are fitted. """ n_photon_dict = dict() for n_photons in range(n_photons_min, n_photons_max, n_photon_step): tau_sim = list() tau_recov = list() n_photons = int(n_photons) for i in range(n_samples): tau = np.random.uniform(tau_min, tau_max) # n_photons = int(5./tau * n_photon_max) param = np.array([tau, gamma, r0, rho]) corrections = np.array([period, g_factor, l1, l2, conv_stop]) model = np.zeros_like(irf_np) bg = np.zeros_like(irf_np) tttrlib.DecayFit23.modelf(param, irf_np, bg, dt, corrections, model) model *= n_photons / np.sum(model) data = np.random.poisson(model) # This performs a with on the data r = fit23(data=data, initial_values=x0, fixed=fixed) # print("tau_sim: %.2f, tau_recov: %s" % (tau, r['x'][0])) tau_sim.append(tau) tau_recov.append(r['x'][0]) n_photon_dict[n_photons] = { 'tau_simulated': np.array(tau_sim), 'tau_recovered': np.array(tau_recov) } devs = list() for k in n_photon_dict: tau_sim = n_photon_dict[k]['tau_simulated'] tau_recov = n_photon_dict[k]['tau_recovered'] dev = (tau_recov - tau_sim) / tau_sim devs.append(dev) """ The figures below demonstrate how well the fluorescence lifetime can be recovered. The left figure displays a typical simulated fluorescence decay. The middle figure illustrated the relative deviation of the recovered fluorescence lifetime :math:`\tau_{recov}` from the simulated fluorescence lifetime :math:`\tau_{sim}`. The figure to the right displays the standard deviation of the relative deviations in dependence of the number of simulated photons. """ fig, ax = plt.subplots(nrows=1, ncols=3, squeeze=True) fig.set_size_inches(12, 4) fig.subplots_adjust(bottom=0.2, left=0.125, right=0.9, wspace=0.3) ax[0].semilogy([x for x in fit23.data], label='Data') ax[0].semilogy([x for x in fit23.irf], label='IRF') ax[0].semilogy([x for x in fit23.model], label='Model') ax[0].set_ylim((0.1, 10000)) ax[0].title.set_text(r'Example decay') ax[0].set_ylabel(r'Counts') ax[0].set_xlabel(r'Channel Nbr.') ax[0].legend() k = list(n_photon_dict.keys())[0] tau_sim = n_photon_dict[k]['tau_simulated'] tau_recov = n_photon_dict[k]['tau_recovered'] dev = (tau_recov - tau_sim) / tau_sim ax[1].title.set_text(r'$\tau_{recov}$ for 20 photons') ax[1].plot(tau_sim, dev, 'o', label='#Photons: %s' % k) ax[1].set_ylim((-1.5, 1.5)) ax[1].set_ylabel(r'$\tau_{recov} - \tau_{sim} / \tau_{sim}$') ax[1].set_xlabel(r'$\tau_{sim}$') sq_dev = np.array(devs)**2 ax[2].plot(list(n_photon_dict.keys()), np.sqrt(sq_dev.mean(axis=1))) ax[2].title.set_text(r'Photon count dependence') ax[2].set_ylabel(r'Std($\tau_{recov} - \tau_{sim} / \tau_{sim}$)') ax[2].set_xlabel(r'Nbr of photons') plt.show() .. rst-class:: sphx-glr-timing **Total running time of the script:** (0 minutes 2.740 seconds) .. _sphx_glr_download_auto_examples_fluorescence_decay_plot_fit23_benchmark_1.py: .. only:: html .. container:: sphx-glr-footer sphx-glr-footer-example .. container:: sphx-glr-download sphx-glr-download-jupyter :download:`Download Jupyter notebook: plot_fit23_benchmark_1.ipynb ` .. container:: sphx-glr-download sphx-glr-download-python :download:`Download Python source code: plot_fit23_benchmark_1.py ` .. only:: html .. rst-class:: sphx-glr-signature `Gallery generated by Sphinx-Gallery `_