from fretbursts import * import os %matplotlib inline import fretbursts.style from IPython.html.widgets import interact, interactive, fixed from IPython.html import widgets from IPython.display import display, display_png, display_svg, clear_output from IPython.core.pylabtools import print_figure import matplotlib as mpl import seaborn as sns sns.set_style('darkgrid') fs = 13 rc={'font.size': fs, 'axes.labelsize': fs, 'legend.fontsize': fs, 'axes.titlesize': fs*1.1, 'xtick.labelsize': fs, 'ytick.labelsize': fs, 'savefig.dpi': 70, #'figure.figsize': (6, 4), #'axes.facecolor': '0.95', 'axes.edgecolor': '0.85', 'grid.color': '0.85', #'axes.linewidth': 1, } sns.set(rc=rc) blue = '#0055d4' green = '#2ca02c' red = "#e74c3c" purple = "#9b59b6" color_brewer = sns.color_palette("Set1", 9) colors = np.array(color_brewer)[(1,0,2,3,4,8,6,7), :] colors = list(colors) colors[:3] = (blue, colors[1], green) sns.set_palette(colors, 8) #sns.palplot(sns.color_palette(colors, 8)) url = 'http://files.figshare.com/1839121/0023uLRpitc_NTP_20dT_0.5GndCl.hdf5' download_file(url, save_dir='./data') file_name = "0023uLRpitc_NTP_20dT_0.5GndCl.hdf5" # Here the folder is the subfolder "data" of current notebook folder folder_name = './data/' full_fname = folder_name + file_name full_fname if os.path.isfile(full_fname): print "Perfect, I found the file!" else: print "Sorry, I can't find the file:\n", full_fname d = loader.hdf5(fname=full_fname) d.add(det_donor_accept=(0, 1), alex_period=4000, D_ON=(2850, 580), A_ON=(900, 2580)) bpl.plot_alternation_hist(d) loader.usalex_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 = Sel(d, select_bursts.size, add_naa=True, th1=30) #dplot(ds, hist_fret) #dplot(ds, hist_S); #ds.E_fitter.fit_histogram(mfit.factory_three_gaussians()) dplot(ds, scatter_alex, figsize=(4,4), mew=1, ms=4, mec='black', color='purple'); dplot(ds, hist2d_alex); def _alex_plot_style(g): """Set plot style and colorbar for an ALEX joint plot. """ g.set_axis_labels(xlabel="E", ylabel="S") g.ax_marg_x.grid(True) g.ax_marg_y.grid(True) plt.setp(g.ax_marg_y.get_xticklabels(), visible=True) plt.setp(g.ax_marg_x.get_yticklabels(), visible=True) g.ax_marg_x.locator_params(axis='y', tight=True, nbins=3) g.ax_marg_y.locator_params(axis='x', tight=True, nbins=3) pos = g.ax_joint.get_position().get_points() X, Y = pos[:, 0], pos[:, 1] cax = plt.axes([1., Y[0], (X[1] - X[0])*0.045, Y[1]-Y[0]]) plt.colorbar(cax=cax) def _hist_bursts(arr, dx, **kwargs): """Wrapper function for calling hist_burst_data() from seaborn plot_marginals(). """ vertical = kwargs.get('vertical', False) data_name = 'S' if vertical else 'E' hist_burst_data(dx, data_name=data_name, **kwargs) def _alex_hexbin_vmax(patches, vmax_fret=True, Smax=0.8): """Return max count E-S hexbin plot in `patches`. When `vmax_fret` is True, returns the max count for S < Smax. Otherwise returns the max count in all the histogram. """ counts = patches.get_array() if vmax_fret: offset = patches.get_offsets() xoffset, yoffset = offset[:, 0], offset[:, 1] mask = yoffset < Smax vmax = counts[mask].max() else: vmax = counts.max() return vmax def _calc_vmin(vmax, vmax_threshold, vmin_default): if vmax <= vmax_threshold: vmin = vmin_default - 0.5*vmax elif vmax_threshold < vmax < 2*vmax_threshold: vmin = vmin_default - 0.5*vmax*((2*vmax_threshold - vmax)/(1.*vmax_threshold)) else: vmin = vmin_default return vmin def alex_jointplot(d, i=0, gridsize=50, cmap='YlGnBu_crop', zero_transparent=True, vmax_fret=True, vmax_threshold=10, vmin_default=0, vmin=None, cmap_compensate=False, joint_kws = {}, marginal_kws = dict(show_kde=True, bandwidth=0.03, binwidth=0.03)): """Plot an ALEX join plot: E-S 2-D histograms and marginal plots. """ g = sns.JointGrid(x=d.E[i], y=d.S[i], ratio=3, space=0.2, xlim=(-0.2, 1.2), ylim=(-0.2, 1.2)) joint_kws.update(gridsize=gridsize, cmap=cmap, extent=(-0.2, 1.2, -0.2, 1.2)) jplot = g.plot_joint(plt.hexbin, mincnt=zero_transparent, **joint_kws) poly = jplot.ax_joint.get_children()[2] vmax = _alex_hexbin_vmax(poly, vmax_fret=vmax_fret) if vmin is None: if cmap_compensate: #vmin = _calc_vmin(vmax, vmax_threshold, vmin_default) vmin = -vmax/3 else: vmin = vmin_default poly.set_clim(vmin, vmax) g.plot_marginals(_hist_bursts, dx=d, **marginal_kws) g.annotate(lambda x, y: x.size, stat='# Bursts', template='{stat}: {val}') _alex_plot_style(g) sns.palplot(sns.color_palette('nipy_spectral', 64)) c = sns.color_palette('nipy_spectral', 64)[2:43] sns.palplot(c) cmap = mpl.colors.LinearSegmentedColormap.from_list('alex', c) cmap.set_under(alpha=0) mpl.cm.register_cmap(name='alex_lv', cmap=cmap) cmap squarehist, _, _ = np.histogram2d(ds.E_, ds.S_, range=((-0.2, 1.2), (-0.2, 1.2)), bins=np.arange(-0.2, 1.2, 0.018)) im = plt.imshow(squarehist[:,::-1].T, extent=(-0.2, 1.2, -0.2, 1.2), vmin=1, cmap='alex_lv', zorder=10) gca().set_axisbelow(True) plt.colorbar() sns.palplot(sns.color_palette('YlGnBu', 64)) c = sns.color_palette('YlGnBu', 64)[16:] sns.palplot(c) cmap = mpl.colors.LinearSegmentedColormap.from_list('YlGnBu_crop', c) cmap.set_under(alpha=0) mpl.cm.register_cmap(name='YlGnBu_crop', cmap=cmap) cmap # Fit E and S to a model and compute KDE bext.bursts_fitter(ds, 'E', binwidth=0.03, bandwidth=0.03, model=mfit.factory_three_gaussians()) bext.bursts_fitter(ds, 'S', binwidth=0.03, bandwidth=0.03, model=mfit.factory_two_gaussians()) sns.set_style('whitegrid') alex_jointplot(ds) alex_jointplot(ds, cmap = 'GnBu', cmap_compensate = True, gridsize = 50, zero_transparent = True, vmax_fret = True, joint_kws = dict(edgecolor='grey')) sns.set_style('darkgrid') alex_jointplot(ds) alex_jointplot(ds, #cmap = 'GnBu', #cmap_compensate = True, gridsize = 50, zero_transparent = True, vmax_fret = True, joint_kws = dict(edgecolor='grey')) sns.set_style('darkgrid') sns.set_style('whitegrid') @interact(cmap_compensate = False, overlay = widgets.RadioButtonsWidget(values=['fit model', 'KDE']), binwidth = widgets.FloatTextWidget(value=0.03, min=0.01, max=1), bandwidth = widgets.FloatTextWidget(value=0.03, min=0.01, max=1), gridsize = (10, 100), min_size=(10, 500, 5), cmap=widgets.DropdownWidget(value='YlGnBu_crop', values=['YlGnBu_crop', 'YlOrRd', 'Blues', 'PuBuGn', 'PuBu', 'GnBu', 'YlGnBu', 'afmhot', 'alex_lv', 'copper', 'summer', 'winter', 'cubehelix']), zero_transparent = True, reverse_cmap = False, vmax_fret = True, ) def plot_(min_size=50, cmap_compensate=False, overlay='KDE', binwidth=0.03, bandwidth=0.03, gridsize=50, cmap='YlGnBu_crop', reverse_cmap=False, vmax_fret=True, zero_transparent=True): dx = Sel(d, select_bursts.size, add_naa=True, th1=min_size) bext.bursts_fitter(dx, 'E', binwidth=binwidth, bandwidth=bandwidth, model=mfit.factory_three_gaussians()) bext.bursts_fitter(dx, 'S', binwidth=binwidth, bandwidth=bandwidth, model=mfit.factory_two_gaussians()) if reverse_cmap: cmap += '_r' if binwidth < 0.01: binwidth = 0.01 if bandwidth < 0.01: bandwidth = 0.01 if overlay == 'fit model': marginal_kws = dict(binwidth=binwidth, show_model=True) else: marginal_kws = dict(binwidth=binwidth, show_kde=True, bandwidth=bandwidth) alex_jointplot(dx, cmap=cmap, gridsize=gridsize, zero_transparent=zero_transparent, vmax_fret=vmax_fret, cmap_compensate=cmap_compensate, marginal_kws=marginal_kws, joint_kws=dict())#edgecolor='grey')) fig = gcf() plt.close() display(fig) from IPython.core.display import HTML HTML(open("./styles/custom2.css", "r").read())