diff --git a/FiCurve.py b/FiCurve.py index 0376809..f9c2249 100644 --- a/FiCurve.py +++ b/FiCurve.py @@ -6,6 +6,7 @@ from scipy.stats import linregress import matplotlib.pyplot as plt from warnings import warn import functions as fu +import helperFunctions as hF class FICurve: @@ -29,12 +30,16 @@ class FICurve: self.f_infinity_fit = [] self.all_calculate_frequency_points() - self.fit_line() - self.fit_boltzmann() + self.f_infinity_fit = hF.fit_clipped_line(self.stimulus_value, self.f_infinities) + self.boltzmann_fit_vars = hF.fit_boltzmann(self.stimulus_value, self.f_zeros) def all_calculate_frequency_points(self): mean_frequencies = self.cell_data.get_mean_isi_frequencies() time_axes = self.cell_data.get_time_axes_mean_frequencies() + stimulus_start = self.cell_data.get_stimulus_start() + stimulus_duration = self.cell_data.get_stimulus_duration() + sampling_interval = self.cell_data.get_sampling_interval() + if len(mean_frequencies) == 0: warn("FICurve:all_calculate_frequency_points(): mean_frequencies is empty.\n" "Was all_calculate_mean_isi_frequencies already called?") @@ -48,110 +53,97 @@ class FICurve: self.f_infinities.append(-1) continue - self.f_zeros.append(self.__calculate_f_zero__(time_axes[i], mean_frequencies[i])) - self.f_baselines.append(self.__calculate_f_baseline__(time_axes[i], mean_frequencies[i])) - self.f_infinities.append(self.__calculate_f_infinity__(time_axes[i], mean_frequencies[i])) - - def fit_line(self): - popt, pcov = curve_fit(fu.clipped_line, self.stimulus_value, self.f_infinities) - self.f_infinity_fit = popt - - def fit_boltzmann(self): - max_f0 = float(max(self.f_zeros)) - min_f0 = 0.1 # float(min(self.f_zeros)) - mean_int = float(np.mean(self.stimulus_value)) - - total_increase = max_f0 - min_f0 - total_change_int = max(self.stimulus_value) - min(self.stimulus_value) - start_k = float((total_increase / total_change_int * 4) / max_f0) - - popt, pcov = curve_fit(fu.full_boltzmann, self.stimulus_value, self.f_zeros, - p0=(max_f0, min_f0, start_k, mean_int), - maxfev=10000, bounds=([0, 0, -np.inf, -np.inf], [5000, 1, np.inf, np.inf])) - - self.boltzmann_fit_vars = popt - - def __calculate_f_baseline__(self, time, frequency, buffer=0.025): - - stim_start = self.cell_data.get_stimulus_start() - time[0] - sampling_interval = self.cell_data.get_sampling_interval() - if stim_start < 0.1: - warn("FICurve:__calculate_f_baseline__(): Quite short delay at the start.") - - start_idx = 0 - end_idx = int((stim_start-buffer)/sampling_interval) - f_baseline = np.mean(frequency[start_idx:end_idx]) - - return f_baseline - - def __calculate_f_zero__(self, time, frequency, peak_buffer_percent=0.05, buffer=0.025): - - stimulus_start = self.cell_data.get_stimulus_start() - time[0] # time start is generally != 0 and != delay - sampling_interval = self.cell_data.get_sampling_interval() - - freq_before = frequency[0:int((stimulus_start - buffer) / sampling_interval)] - min_before = min(freq_before) - max_before = max(freq_before) - mean_before = np.mean(freq_before) - - # time where the f-zero is searched in - start_idx = int((stimulus_start-0.1*buffer) / sampling_interval) - end_idx = int((stimulus_start + buffer) / sampling_interval) - - min_during_start_of_stim = min(frequency[start_idx:end_idx]) - max_during_start_of_stim = max(frequency[start_idx:end_idx]) - - if abs(mean_before-min_during_start_of_stim) > abs(max_during_start_of_stim-mean_before): - f_zero = min_during_start_of_stim - else: - f_zero = max_during_start_of_stim - - peak_buffer = (max_before - min_before) * peak_buffer_percent - if min_before - peak_buffer <= f_zero <= max_before + peak_buffer: - end_idx = start_idx + int((end_idx-start_idx)/2) - f_zero = np.mean(frequency[start_idx:end_idx]) - - return f_zero - - # start_idx = int(stimulus_start / sampling_interval) - # end_idx = int((stimulus_start + buffer*2) / sampling_interval) - # - # freq_before = frequency[start_idx-(int(length_of_mean/sampling_interval)):start_idx] - # fb_mean = np.mean(freq_before) - # fb_std = np.std(freq_before) - # - # peak_frequency = fb_mean - # count = 0 - # for i in range(start_idx + 1, end_idx): - # if fb_mean-3*fb_std <= frequency[i] <= fb_mean+3*fb_std: - # continue - # - # if abs(frequency[i] - fb_mean) > abs(peak_frequency - fb_mean): - # peak_frequency = frequency[i] - # count += 1 - - # return peak_frequency - - def __calculate_f_infinity__(self, time, frequency, length=0.1, buffer=0.025): - stimulus_end_time = self.cell_data.get_stimulus_start() + self.cell_data.get_stimulus_duration() - time[0] - - start_idx = int((stimulus_end_time - length - buffer) / self.cell_data.get_sampling_interval()) - end_idx = int((stimulus_end_time - buffer) / self.cell_data.get_sampling_interval()) - - # TODO add way to plot detected f_zero, f_inf, f_base. With detection of remaining slope? - # x = np.arange(start_idx, end_idx, 1) # time[start_idx:end_idx] - # slope, intercept, r_value, p_value, std_err = linregress(x, frequency[start_idx:end_idx]) - # if p_value < 0.0001: - # plt.title("significant slope: {:.2f}, p: {:.5f}, r: {:.5f}".format(slope, p_value, r_value)) - # plt.plot(x, [i*slope + intercept for i in x], color="black") - # - # - # plt.plot((start_idx, end_idx), (np.mean(frequency[start_idx:end_idx]), np.mean(frequency[start_idx:end_idx])), label="f_inf") - # plt.legend() - # plt.show() - # plt.close() - - return np.mean(frequency[start_idx:end_idx]) + f_zero = hF.detect_f_zero_in_frequency_trace(time_axes[i], mean_frequencies[i], + stimulus_start, sampling_interval) + self.f_zeros.append(f_zero) + f_baseline = hF.detect_f_baseline_in_freq_trace(time_axes[i], mean_frequencies[i], + stimulus_start, sampling_interval) + self.f_baselines.append(f_baseline) + f_infinity = hF.detect_f_infinity_in_freq_trace(time_axes[i], mean_frequencies[i], + stimulus_start, stimulus_duration, sampling_interval) + self.f_infinities.append(f_infinity) + + # def __calculate_f_baseline__(self, time, frequency, buffer=0.025): + # + # stim_start = self.cell_data.get_stimulus_start() - time[0] + # sampling_interval = self.cell_data.get_sampling_interval() + # if stim_start < 0.1: + # warn("FICurve:__calculate_f_baseline__(): Quite short delay at the start.") + # + # start_idx = 0 + # end_idx = int((stim_start-buffer)/sampling_interval) + # f_baseline = np.mean(frequency[start_idx:end_idx]) + # + # return f_baseline + # + # def __calculate_f_zero__(self, time, frequency, peak_buffer_percent=0.05, buffer=0.025): + # + # stimulus_start = self.cell_data.get_stimulus_start() - time[0] # time start is generally != 0 and != delay + # sampling_interval = self.cell_data.get_sampling_interval() + # + # freq_before = frequency[0:int((stimulus_start - buffer) / sampling_interval)] + # min_before = min(freq_before) + # max_before = max(freq_before) + # mean_before = np.mean(freq_before) + # + # # time where the f-zero is searched in + # start_idx = int((stimulus_start-0.1*buffer) / sampling_interval) + # end_idx = int((stimulus_start + buffer) / sampling_interval) + # + # min_during_start_of_stim = min(frequency[start_idx:end_idx]) + # max_during_start_of_stim = max(frequency[start_idx:end_idx]) + # + # if abs(mean_before-min_during_start_of_stim) > abs(max_during_start_of_stim-mean_before): + # f_zero = min_during_start_of_stim + # else: + # f_zero = max_during_start_of_stim + # + # peak_buffer = (max_before - min_before) * peak_buffer_percent + # if min_before - peak_buffer <= f_zero <= max_before + peak_buffer: + # end_idx = start_idx + int((end_idx-start_idx)/2) + # f_zero = np.mean(frequency[start_idx:end_idx]) + # + # return f_zero + # + # # start_idx = int(stimulus_start / sampling_interval) + # # end_idx = int((stimulus_start + buffer*2) / sampling_interval) + # # + # # freq_before = frequency[start_idx-(int(length_of_mean/sampling_interval)):start_idx] + # # fb_mean = np.mean(freq_before) + # # fb_std = np.std(freq_before) + # # + # # peak_frequency = fb_mean + # # count = 0 + # # for i in range(start_idx + 1, end_idx): + # # if fb_mean-3*fb_std <= frequency[i] <= fb_mean+3*fb_std: + # # continue + # # + # # if abs(frequency[i] - fb_mean) > abs(peak_frequency - fb_mean): + # # peak_frequency = frequency[i] + # # count += 1 + # + # # return peak_frequency + # + # def __calculate_f_infinity__(self, time, frequency, length=0.1, buffer=0.025): + # stimulus_end_time = self.cell_data.get_stimulus_start() + self.cell_data.get_stimulus_duration() - time[0] + # + # start_idx = int((stimulus_end_time - length - buffer) / self.cell_data.get_sampling_interval()) + # end_idx = int((stimulus_end_time - buffer) / self.cell_data.get_sampling_interval()) + # + # # TODO add way to plot detected f_zero, f_inf, f_base. With detection of remaining slope? + # # x = np.arange(start_idx, end_idx, 1) # time[start_idx:end_idx] + # # slope, intercept, r_value, p_value, std_err = linregress(x, frequency[start_idx:end_idx]) + # # if p_value < 0.0001: + # # plt.title("significant slope: {:.2f}, p: {:.5f}, r: {:.5f}".format(slope, p_value, r_value)) + # # plt.plot(x, [i*slope + intercept for i in x], color="black") + # # + # # + # # plt.plot((start_idx, end_idx), (np.mean(frequency[start_idx:end_idx]), np.mean(frequency[start_idx:end_idx])), label="f_inf") + # # plt.legend() + # # plt.show() + # # plt.close() + # + # return np.mean(frequency[start_idx:end_idx]) def get_f_zero_inverse_at_frequency(self, frequency): b_vars = self.boltzmann_fit_vars @@ -226,7 +218,6 @@ class FICurve: if comp_f_zeros is not None: plt.plot(self.stimulus_value, comp_f_zeros, 'o', color='wheat', label='comp_values f_zero') - plt.legend() plt.ylabel("Frequency [Hz]") if self.using_contrast: @@ -236,6 +227,7 @@ class FICurve: if savepath is None: plt.show() else: + print("save") plt.savefig(savepath + "fi_curve.png") plt.close() @@ -246,9 +238,9 @@ class FICurve: for i in range(len(mean_frequencies)): fig, axes = plt.subplots(1, 1, sharex="all") axes.plot(time_axes[i], mean_frequencies[i], label="voltage") - axes.plot((time_axes[i][0],time_axes[i][-1]), (self.f_zeros[i], self.f_zeros[i]), label="f_zero") - axes.plot((time_axes[i][0],time_axes[i][-1]), (self.f_infinities[i], self.f_infinities[i]), '--', label="f_inf") - axes.plot((time_axes[i][0],time_axes[i][-1]), (self.f_baselines[i], self.f_baselines[i]), label="f_base") + axes.plot((time_axes[i][0], time_axes[i][-1]), (self.f_zeros[i], self.f_zeros[i]), label="f_zero") + axes.plot((time_axes[i][0], time_axes[i][-1]), (self.f_infinities[i], self.f_infinities[i]), '--', label="f_inf") + axes.plot((time_axes[i][0], time_axes[i][-1]), (self.f_baselines[i], self.f_baselines[i]), label="f_base") axes.set_title(str(self.stimulus_value[i])) plt.legend()