diff --git a/fit_lifacnoise.py b/fit_lifacnoise.py index 59bd5c3..661f10b 100644 --- a/fit_lifacnoise.py +++ b/fit_lifacnoise.py @@ -126,7 +126,7 @@ class Fitter: # only eod with amplitude 1 and no modulation _, spiketimes = self.model.simulate_fast(base_stimulus, 30) - baseline_freq = hF.mean_freq_of_spiketimes_after_time_x(spiketimes, freq_sampling_rate, 5) + baseline_freq = hF.mean_freq_of_spiketimes_after_time_x(spiketimes, 5) # print("model:", baseline_freq, "data:", self.baseline_freq) relative_spiketimes = np.array([s % (1/self.eod_freq) for s in spiketimes if s > 0]) @@ -142,7 +142,7 @@ class Fitter: if len(spiketimes) < 2: f_infinities.append(0) else: - f_infinity = hF.mean_freq_of_spiketimes_after_time_x(spiketimes, freq_sampling_rate, 0.5) + f_infinity = hF.mean_freq_of_spiketimes_after_time_x(spiketimes, 0.5) f_infinities.append(f_infinity) popt, pcov = curve_fit(fu.line, self.fi_contrasts, f_infinities, maxfev=10000) diff --git a/helperFunctions.py b/helperFunctions.py index 74a0d2f..5a35ba6 100644 --- a/helperFunctions.py +++ b/helperFunctions.py @@ -155,19 +155,33 @@ def calculate_mean_of_frequency_traces(trial_time_traces, trial_frequency_traces return shortened_time, mean_freq -def mean_freq_of_spiketimes_after_time_x(spiketimes, sampling_interval, time_x, time_in_ms=False): +def mean_freq_of_spiketimes_after_time_x(spiketimes, time_x, time_in_ms=False): """ Calculates the mean frequency of the portion of spiketimes that is after last_x_time """ if len(spiketimes) <= 1: return 0 - freq = calculate_isi_frequency_trace(spiketimes, sampling_interval, time_in_ms) - # returned frequency starts at the - idx = int((time_x-spiketimes[0]) / sampling_interval) - rest_array = freq[idx:] - mean_freq = np.mean(rest_array) + relevant_spikes = spiketimes[spiketimes > time_x] + + if len(relevant_spikes) <= 1: + return 0 + if time_in_ms: + relevant_spikes = relevant_spikes / 1000 + isis = np.diff(relevant_spikes) + isi_freqs = 1 / isis + weights = isis / min(isis) + + mean_freq = sum(isi_freqs * weights) / sum(weights) + return mean_freq + # freq = calculate_isi_frequency_trace(spiketimes, sampling_interval, time_in_ms) + # # returned frequency starts at the + # idx = int((time_x-spiketimes[0]) / sampling_interval) + # rest_array = freq[idx:] + # mean_freq = np.mean(rest_array) + # return mean_freq + def calculate_mean_isi_freq(spiketimes, time_in_ms=False): if len(spiketimes) < 2: diff --git a/models/LIFACnoise.py b/models/LIFACnoise.py index 40b3a11..8ba3618 100644 --- a/models/LIFACnoise.py +++ b/models/LIFACnoise.py @@ -163,7 +163,7 @@ class LifacNoiseModel(AbstractModel): base_stimulus = SinusAmplitudeModulationStimulus(base_stimulus_freq, 0, 0) _, spiketimes = self.simulate_fast(base_stimulus, 30) time_x = 5 - baseline_freq = hF.mean_freq_of_spiketimes_after_time_x(spiketimes, self.get_sampling_interval(), time_x) + baseline_freq = hF.mean_freq_of_spiketimes_after_time_x(spiketimes, time_x) relative_spiketimes = np.array([s % (1 / base_stimulus_freq) for s in spiketimes]) eod_durations = np.full((len(spiketimes)), 1 / base_stimulus_freq) @@ -184,7 +184,7 @@ class LifacNoiseModel(AbstractModel): stimulus = SinusAmplitudeModulationStimulus(base_freq, contrast, modulation_frequency) _, spiketimes = self.simulate_fast(stimulus, 1) - f_infinity = hF.mean_freq_of_spiketimes_after_time_x(spiketimes, self.get_sampling_interval(), 0.3) + f_infinity = hF.mean_freq_of_spiketimes_after_time_x(spiketimes, 0.3) f_infinities.append(f_infinity) popt, pcov = curve_fit(fu.line, contrasts, f_infinities, maxfev=10000) @@ -239,7 +239,7 @@ def test_v_offset(model: LifacNoiseModel, v_offset, base_stimulus, simulation_le model.set_variable("v_offset", v_offset) _, spiketimes = model.simulate_fast(base_stimulus, simulation_length) - freq = hF.mean_freq_of_spiketimes_after_time_x(spiketimes, 0.0005, simulation_length/3) + freq = hF.mean_freq_of_spiketimes_after_time_x(spiketimes, simulation_length / 3) return freq diff --git a/unittests/testFrequencyFunctions.py b/unittests/testFrequencyFunctions.py index ed342ce..bb03e45 100644 --- a/unittests/testFrequencyFunctions.py +++ b/unittests/testFrequencyFunctions.py @@ -31,7 +31,7 @@ class FrequencyFunctionsTester(unittest.TestCase): for freq in self.frequencies: for n in self.noise_levels: spikes = generate_jittered_spiketimes(freq, n, end=simulation_time) - sim_freq = hF.mean_freq_of_spiketimes_after_time_x(spikes, 0.00005, simulation_time/4, time_in_ms=False) + sim_freq = hF.mean_freq_of_spiketimes_after_time_x(spikes, simulation_time / 4, time_in_ms=False) max_diff = round(n*(10+0.7*np.sqrt(freq)), 2) # print("noise: {:.2f}".format(n), "\texpected: {:.2f}".format(freq), "\tgotten: {:.2f}".format(round(sim_freq, 2)), "\tfreq diff: {:.2f}".format(abs(freq-round(sim_freq, 2))), "\tmax_diff:", max_diff) @@ -223,7 +223,7 @@ def test_distribution(): print("#### - freq:", freq, "noise level:", n ) for reps in range(repetitions): spikes = generate_jittered_spiketimes(freq, n, end=simulation_time) - sim_freq = hF.mean_freq_of_spiketimes_after_time_x(spikes, 0.0002, simulation_time / 4, time_in_ms=False) + sim_freq = hF.mean_freq_of_spiketimes_after_time_x(spikes, simulation_time / 4, time_in_ms=False) diffs.append(sim_freq-freq) diffs_per_noise.append(diffs) diff --git a/unittests/testLifacNoise.py b/unittests/testLifacNoise.py index e78d4c9..3b96afd 100644 --- a/unittests/testLifacNoise.py +++ b/unittests/testLifacNoise.py @@ -55,7 +55,7 @@ class HelperFunctionsTester(unittest.TestCase): test_model.set_variable("v_offset", offset) _, spikes = test_model.simulate_fast(stimulus, 5) - goal_freq = hF.mean_freq_of_spiketimes_after_time_x(spikes, 0.0005, 1) + goal_freq = hF.mean_freq_of_spiketimes_after_time_x(spikes, 1) if goal_freq <= threshold: print("test Offset ({:.1f}) generates a too low frequency: {:.2f}".format(offset, goal_freq))