diff --git a/Figures_Baseline.py b/Figures_Baseline.py index 6357f1f..e193fe7 100644 --- a/Figures_Baseline.py +++ b/Figures_Baseline.py @@ -9,12 +9,9 @@ import helperFunctions as hF from CellData import CellData from Baseline import BaselineCellData -from FiCurve import FICurveCellData +from FiCurve import FICurveCellData, FICurveModel import Figure_constants as consts -MODEL_COLOR = "orange" -DATA_COLOR = "blue" - -DATA_SAVE_PATH = "data/figure_data/" +from ModelFit import get_best_fit EXAMPLE_CELL = "data/final/2012-12-20-ac-invivo-1" @@ -25,13 +22,68 @@ def main(): # data_fi_curve() p_unit_example() fi_point_detection() + + + # test_fi_curve_colors() pass -def p_unit_example(): +def test_fi_curve_colors(): + example_cell_fit = "results/final_2/2012-12-20-ac-invivo-1" cell = EXAMPLE_CELL + cell_data = CellData(cell) + fit = get_best_fit(example_cell_fit) + + fig, axes = plt.subplots(1, 3) + + axes[0].set_title("Cell") + fi_curve = FICurveCellData(cell_data, cell_data.get_fi_contrasts(), save_dir=cell_data.get_data_path()) + contrasts = cell_data.get_fi_contrasts() + f_zeros = fi_curve.get_f_zero_frequencies() + f_infs = fi_curve.get_f_inf_frequencies() + axes[0].plot(contrasts, f_zeros, ',', marker=consts.f0_marker, color=consts.COLOR_DATA_f0) + axes[0].plot(contrasts, f_infs, ',', marker=consts.finf_marker, color=consts.COLOR_DATA_finf) + + x_values = np.arange(min(contrasts), max(contrasts), (max(contrasts) - min(contrasts)) / 1000) + f_inf_fit = fi_curve.f_inf_fit + f_zero_fit = fi_curve.f_zero_fit + f_zero_fit = [fu.full_boltzmann(x, f_zero_fit[0], f_zero_fit[1], f_zero_fit[2], f_zero_fit[3]) for x in x_values] + f_inf_fit = [fu.clipped_line(x, f_inf_fit[0], f_inf_fit[1]) for x in x_values] + axes[0].plot(x_values, f_zero_fit, color=consts.COLOR_DATA_f0) + axes[0].plot(x_values, f_inf_fit, color=consts.COLOR_DATA_finf) + + axes[2].plot(x_values, f_zero_fit, color=consts.COLOR_DATA_f0) + axes[2].plot(x_values, f_inf_fit, color=consts.COLOR_DATA_finf) + + axes[1].set_title("Model") + model = fit.get_model() + fi_curve = FICurveModel(model, contrasts, eod_frequency=cell_data.get_eod_frequency()) + + f_zeros = fi_curve.get_f_zero_frequencies() + f_infs = fi_curve.get_f_inf_frequencies() + axes[1].plot(contrasts, f_zeros, ',', marker=consts.f0_marker, color=consts.COLOR_MODEL_f0) + axes[1].plot(contrasts, f_infs, ',', marker=consts.finf_marker, color=consts.COLOR_MODEL_finf) + + x_values = np.arange(min(contrasts), max(contrasts), (max(contrasts) - min(contrasts)) / 1000) + f_inf_fit = fi_curve.f_inf_fit + f_zero_fit = fi_curve.f_zero_fit + f_zero_fit = [fu.full_boltzmann(x, f_zero_fit[0], f_zero_fit[1], f_zero_fit[2], f_zero_fit[3]) for x in x_values] + f_inf_fit = [fu.clipped_line(x, f_inf_fit[0], f_inf_fit[1]) for x in x_values] + axes[1].plot(x_values, f_zero_fit, color=consts.COLOR_MODEL_f0) + axes[1].plot(x_values, f_inf_fit, color=consts.COLOR_MODEL_finf) + + axes[2].plot(contrasts, f_zeros, ",", marker=consts.f0_marker, color=consts.COLOR_MODEL_f0) + axes[2].plot(contrasts, f_infs, ",", marker=consts.finf_marker, color=consts.COLOR_MODEL_finf) + + plt.show() + plt.close() + +def p_unit_example(): + cell = EXAMPLE_CELL + cell_data = CellData(cell) + print("p-unit example eodf:", cell_data.get_eod_frequency()) base = BaselineCellData(cell_data) base.load_values(cell_data.get_data_path()) print("burstiness of example cell:", base.get_burstiness()) @@ -51,12 +103,13 @@ def p_unit_example(): start = 0 duration = 0.10 - ax.plot(np.array(time[:int(duration/step)]) - start, v1[:int(duration/step)]) - ax.eventplot([s for s in spiketimes if start < s < start + duration], lineoffsets=max(v1[:int(duration/step)])+0.75, color="black") + ax.plot((np.array(time[:int(duration/step)]) - start) * 1000, v1[:int(duration/step)], consts.COLOR_DATA) + ax.eventplot([s * 1000 for s in spiketimes if start < s < start + duration], lineoffsets=max(v1[:int(duration/step)])+1.25, + color="black", linelengths=2) ax.set_ylabel('Voltage in mV') - ax.set_xlabel('Time in s') + ax.set_xlabel('Time in ms') ax.set_title("Baseline Firing") - ax.set_xlim((0, duration)) + ax.set_xlim((0, duration*1000)) # ISI-hist @@ -66,19 +119,19 @@ def p_unit_example(): isi = np.array(base.get_interspike_intervals()) / eod_period # ISI in ms maximum = max(isi) bins = np.arange(0, maximum * 1.01, 0.1) - ax.hist(isi, bins=bins) - ax.set_ylabel("Count") + ax.hist(isi, bins=bins, color=consts.COLOR_DATA, density=True) + ax.set_ylabel("Density") ax.set_xlabel("ISI in EOD periods") ax.set_title("ISI Histogram") # Serial correlation - ax = fig.add_subplot(gs[2, 0]) + ax = fig.add_subplot(gs[1, 1]) sc = base.get_serial_correlation(10) ax.plot(range(11), [0 for _ in range(11)], color="darkgrey", alpha=0.8) - ax.plot(range(11), [1] + list(sc)) - ax.plot(range(11), [1] + list(sc), '+') + ax.plot(range(11), [1] + list(sc), color=consts.COLOR_DATA) + ax.plot(range(11), [1] + list(sc), '+', color="black") ax.set_xlabel("Lag") ax.set_ylabel("SC") @@ -90,14 +143,14 @@ def p_unit_example(): # FI-Curve trace - ax = fig.add_subplot(gs[1, 1]) + ax = fig.add_subplot(gs[2, 0]) f_trace_times, f_traces = fi.get_mean_time_and_freq_traces() part = 0.4 + 0.2 + 0.2 # stim duration + delay up front and a part of the "delay" at the back idx = int(part/step) - ax.plot(f_trace_times[-1][:idx], f_traces[-1][:idx]) + ax.plot(f_trace_times[-1][:idx], f_traces[-1][:idx], color=consts.COLOR_DATA) # strength = 200 # smoothed = np.convolve(f_traces[-1][:idx], np.ones(strength)/strength) # ax.plot(f_trace_times[-1][:idx], smoothed[int(strength/2):idx + int(strength/2)]) @@ -116,22 +169,28 @@ def p_unit_example(): f_zeros = fi.get_f_zero_frequencies() f_infties = fi.get_f_inf_frequencies() - ax.plot(contrasts, f_zeros, 'x') - ax.plot(contrasts, f_infties, 'o') + ax.plot(contrasts, f_zeros, ',', marker=consts.f0_marker, color=consts.COLOR_DATA_f0) + ax.plot(contrasts, f_infties, ',', marker=consts.finf_marker, color=consts.COLOR_DATA_finf) x_values = np.arange(min(contrasts), max(contrasts) + 0.0001, (max(contrasts)-min(contrasts)) / 1000) f_zero_fit = [fu.full_boltzmann(x, fi.f_zero_fit[0], fi.f_zero_fit[1], fi.f_zero_fit[2], fi.f_zero_fit[3]) for x in x_values] f_inf_fit = [fu.clipped_line(x, fi.f_inf_fit[0], fi.f_inf_fit[1]) for x in x_values] - ax.plot(x_values, f_zero_fit) - ax.plot(x_values, f_inf_fit) + ax.plot(x_values, f_zero_fit, color=consts.COLOR_DATA_f0) + ax.plot(x_values, f_inf_fit, color=consts.COLOR_DATA_finf) # ax.set_xlim((0, 10)) # ax.set_ylim((-1, 1)) ax.set_xlabel("Contrast") ax.set_ylabel("Frequency in Hz") + ax.set_xticks([-0.2, -0.1, 0, 0.1, 0.2]) + ax.set_xlim((-0.21, 0.2)) + ylim = ax.get_ylim() + ax.set_ylim((0, ylim[1])) ax.set_title("f-I Curve") plt.tight_layout() + consts.set_figure_labels(xoffset=-2.5, yoffset=2.2) + fig.label_axes() plt.savefig("thesis/figures/p_unit_example.png", transparent=True) plt.close() @@ -154,10 +213,13 @@ def fi_point_detection(): f_inf = fi.get_f_inf_frequencies()[-1] f_inf_idx = fi.indices_f_inf[-1] + f_baseline = fi.get_f_baseline_frequencies()[-1] + f_base_idx = fi.indices_f_baseline[-1] axes[0].plot(f_trace_times[-1][:idx], f_traces[-1][:idx]) - axes[0].plot([f_trace_times[-1][idx] for idx in f_zero_idx], (f_zero, ), "o") - axes[0].plot([f_trace_times[-1][idx] for idx in f_inf_idx], (f_inf, f_inf), color="orange", linewidth=4) + axes[0].plot([f_trace_times[-1][idx] for idx in f_zero_idx], (f_zero, ), ",", marker=consts.f0_marker, color=consts.COLOR_DATA_f0) + axes[0].plot([f_trace_times[-1][idx] for idx in f_inf_idx], (f_inf, f_inf), color=consts.COLOR_DATA_finf, linewidth=4) + axes[0].plot([f_trace_times[-1][idx] for idx in f_base_idx], (f_baseline, f_baseline), color="grey", linewidth=4) # mark stim start and end: stim_start = cell_data.get_stimulus_start() @@ -176,15 +238,15 @@ def fi_point_detection(): f_zeros = fi.get_f_zero_frequencies() f_infties = fi.get_f_inf_frequencies() - axes[1].plot(contrasts, f_zeros, 'x') - axes[1].plot(contrasts, f_infties, 'o') + axes[1].plot(contrasts, f_zeros, ",", marker=consts.f0_marker, color=consts.COLOR_DATA_f0) + axes[1].plot(contrasts, f_infties, ",", marker=consts.finf_marker, color=consts.COLOR_DATA_finf) x_values = np.arange(min(contrasts), max(contrasts) + 0.0001, (max(contrasts) - min(contrasts)) / 1000) f_zero_fit = [fu.full_boltzmann(x, fi.f_zero_fit[0], fi.f_zero_fit[1], fi.f_zero_fit[2], fi.f_zero_fit[3]) for x in x_values] f_inf_fit = [fu.clipped_line(x, fi.f_inf_fit[0], fi.f_inf_fit[1]) for x in x_values] - axes[1].plot(x_values, f_zero_fit) - axes[1].plot(x_values, f_inf_fit) + axes[1].plot(x_values, f_zero_fit, color=consts.COLOR_DATA_f0) + axes[1].plot(x_values, f_inf_fit, color=consts.COLOR_DATA_finf) axes[1].set_xlabel("Contrast in %") # axes[1].set_ylabel("Frequency in Hz") @@ -193,6 +255,8 @@ def fi_point_detection(): plt.tight_layout() + consts.set_figure_labels(xoffset=-2.5) + fig.label_axes() plt.savefig("thesis/figures/f_point_detection.png", transparent=True) plt.close() diff --git a/Figures_results.py b/Figures_results.py index 4922395..37bfd40 100644 --- a/Figures_results.py +++ b/Figures_results.py @@ -1,18 +1,56 @@ import numpy as np import matplotlib.pyplot as plt from analysis import get_fit_info, get_behaviour_values, calculate_percent_errors +from ModelFit import get_best_fit +from Baseline import BaselineModel, BaselineCellData +import Figure_constants as consts def main(): - dir_path = "results/final_1/" + dir_path = "results/final_2/" fits_info = get_fit_info(dir_path) - cell_behaviour, model_behaviour = get_behaviour_values(fits_info) + # cell_behaviour, model_behaviour = get_behaviour_values(fits_info) # behaviour_overview_pairs(cell_behaviour, model_behaviour) # errors = calculate_percent_errors(fits_info) # create_boxplots(errors) + example_good_hist_fits(dir_path) + + +def example_good_hist_fits(dir_path): + strong_bursty_cell = "2018-05-08-ac-invivo-1" + bursty_cell = "2014-03-19-ad-invivo-1" + non_bursty_cell = "2012-12-21-am-invivo-1" + + fig, axes = plt.subplots(1, 3, sharex="all", figsize=consts.FIG_SIZE_MEDIUM_WIDE) + + for i, cell in enumerate([non_bursty_cell, bursty_cell, strong_bursty_cell]): + fit_dir = dir_path + cell + "/" + fit = get_best_fit(fit_dir) + + cell_data = fit.get_cell_data() + eodf = cell_data.get_eod_frequency() + + model = fit.get_model() + baseline_model = BaselineModel(model, eodf, trials=5) + + model_isi = np.array(baseline_model.get_interspike_intervals()) * eodf + cell_isi = BaselineCellData(cell_data).get_interspike_intervals() * eodf + + bins = np.arange(0, 0.025, 0.0001) * eodf + axes[i].hist(cell_isi, bins=bins, density=True, alpha=0.5, color=consts.COLOR_DATA) + axes[i].hist(model_isi, bins=bins, density=True, alpha=0.5, color=consts.COLOR_MODEL) + axes[i].set_xlabel("ISI in EOD periods") + axes[0].set_ylabel("Density") + plt.tight_layout() + consts.set_figure_labels(xoffset=-2.5) + fig.label_axes() + + plt.savefig(consts.SAVE_FOLDER + "example_good_isi_hist_fits.png", transparent=True) + plt.close() + def create_boxplots(errors): @@ -58,11 +96,21 @@ def overview_pair(cell, model, titles): wspace=0.2, hspace=0.05) for i in range(len(cell)): + if titles[i] == "f_zero_slope": + length_before = len(cell[i]) + idx = np.array(cell[i]) < 30000 + cell[i] = np.array(cell[i])[idx] + model[i] = np.array(model[i])[idx] + + idx = np.array(model[i]) < 30000 + cell[i] = np.array(cell[i])[idx] + model[i] = np.array(model[i])[idx] + print("removed {} values from f_zero_slope plot.".format(length_before - len(cell[i]))) ax = fig.add_subplot(gs[1, i]) ax_histx = fig.add_subplot(gs[0, i], sharex=ax) scatter_hist(cell[i], model[i], ax, ax_histx, titles[i]) - plt.tight_layout() + # plt.tight_layout() plt.show() @@ -88,7 +136,7 @@ def grouped_error_overview_behaviour_dist(cell_behaviours, model_behaviours): # use the previously defined function scatter_hist(cell_behaviours[behaviour], model_behaviours[behaviour], ax, ax_histx, behaviour) - # plt.tight_layout() + plt.tight_layout() plt.show() @@ -96,9 +144,8 @@ def scatter_hist(cell_values, model_values, ax, ax_histx, behaviour, ax_histy=No # copied from matplotlib # no labels - ax_histx.tick_params(axis="cell_values", labelbottom=False) + ax_histx.tick_params(axis="cell", labelbottom=False) # ax_histy.tick_params(axis="model_values", labelleft=False) - # the scatter plot: ax.scatter(cell_values, model_values) @@ -109,17 +156,13 @@ def scatter_hist(cell_values, model_values, ax, ax_histx, behaviour, ax_histy=No ax.set_xlabel("cell") ax.set_ylabel("model") - # now determine nice limits by hand: - binwidth = 0.25 - xymax = max(np.max(np.abs(cell_values)), np.max(np.abs(model_values))) - lim = (int(xymax/binwidth) + 1) * binwidth - - bins = np.arange(-lim, lim + binwidth, binwidth) ax_histx.hist(cell_values, color="blue", alpha=0.5) ax_histx.hist(model_values, color="orange", alpha=0.5) - # ax_histx.set_xticklabels([]) - # ax_histx.set_xticks(ax.get_xticks()) - # ax_histx.set_xlim(ax.get_xlim()) + ax_labels = ax.get_xticklabels() + ax_histx.set_xticklabels([]) + ax.set_xticklabels(ax_labels) + ax_histx.set_xticks(ax.get_xticks()) + ax_histx.set_xlim(ax.get_xlim()) ax_histx.set_title(behaviour) # ax_histy.hist(y, bins=bins, orientation='horizontal') diff --git a/thesis/figures/example_good_isi_hist_fits.png b/thesis/figures/example_good_isi_hist_fits.png new file mode 100644 index 0000000..15a5e08 Binary files /dev/null and b/thesis/figures/example_good_isi_hist_fits.png differ diff --git a/thesis/figures/f_point_detection.png b/thesis/figures/f_point_detection.png index d05aea3..f61cb50 100644 Binary files a/thesis/figures/f_point_detection.png and b/thesis/figures/f_point_detection.png differ diff --git a/thesis/figures/p_unit_example.png b/thesis/figures/p_unit_example.png index 98e34ad..3dac30a 100644 Binary files a/thesis/figures/p_unit_example.png and b/thesis/figures/p_unit_example.png differ