diff --git a/programming/exercises/STA/p-unit_spike_times.mat b/programming/exercises/STA/p-unit_spike_times.mat new file mode 100644 index 0000000..75f53e6 Binary files /dev/null and b/programming/exercises/STA/p-unit_spike_times.mat differ diff --git a/programming/exercises/STA/p-unit_stimulus.mat b/programming/exercises/STA/p-unit_stimulus.mat new file mode 100644 index 0000000..1e88d5e Binary files /dev/null and b/programming/exercises/STA/p-unit_stimulus.mat differ diff --git a/programming/exercises/STA/sta.m b/programming/exercises/STA/sta.m new file mode 100644 index 0000000..50df6ff --- /dev/null +++ b/programming/exercises/STA/sta.m @@ -0,0 +1,18 @@ +function [sta, std_sta, valid_spikes]= sta(stimulus, spike_times, count, sampling_rate) + +snippets = zeros(numel(spike_times), 2*count); +valid_spikes = 1; +for i = 1:numel(spike_times) + t = spike_times(i); + index = round(t*sampling_rate); + if index < count || (index + count) > length(stimulus) + continue + end + snippets(valid_spikes,:) = stimulus(index-count:index+count-1); + valid_spikes = valid_spikes + 1; +end + +snippets(end-(end-valid_spikes):end,:) = []; + +sta = mean(snippets, 1); +std_sta = std(snippets,[],1); \ No newline at end of file diff --git a/programming/exercises/STA/sta_script.m b/programming/exercises/STA/sta_script.m new file mode 100644 index 0000000..a5b6799 --- /dev/null +++ b/programming/exercises/STA/sta_script.m @@ -0,0 +1,62 @@ +clear +clc + +%% load data +load p-unit_spike_times.mat +load p-unit_stimulus.mat +sample_rate = 20000; + +%% calculate STA +all_times = []; +for i = 1:length(spike_times) + all_times = cat(1, all_times, spike_times{i}); +end + +[sta, sta_sd, num] = sta(stimulus_strong(:,2), all_times, 1000, sample_rate); +fig = figure(); +set(fig, 'PaperUnits', 'centimeters'); +set(fig, 'PaperSize', [11.7 9.0]); +set(fig, 'PaperPosition',[0.0 0.0 11.7 9.0]); +set(fig,'Color', 'white') +plot(((1:length(sta))-1000)/sample_rate, sta) +xlabel('time [s]') +ylabel('stimulus') + +%% reverse reconstruction + +% make binary representation of the spike times +binary_spikes = zeros(size(stimulus_strong, 1), length(spike_times)); +estimated_stims = zeros(size(binary_spikes)); +for i = 1:length(spike_times) + binary_spikes(round(spike_times{i}*sample_rate), i) = 1; + estimated_stims(:,i) = conv(binary_spikes(:,i), sta, 'same'); +end + +fig = figure(); +set(fig, 'PaperUnits', 'centimeters'); +set(fig, 'PaperSize', [11.7 9.0]); +set(fig, 'PaperPosition',[0.0 0.0 11.7 9.0]); +set(fig,'Color', 'white') +plot(stimulus_strong(:,1), stimulus_strong(:,2), 'displayname','original') +hold on +plot(stimulus_strong(:,1), mean(estimated_stims,2), 'r', 'displayname', 'reconstruction') +xlabel('time [s]') +ylabel('stimulus') +legend show + +%% calculate STC + +% we need to downsample the data otherwise the covariance matrixs gets too +% large 20Khz to 1kHz + +downsampled_binary = zeros(size(stimulus_strong, 1)/20, length(spike_times)); +downsampled_stim = zeros(size(downsampled_binary,1),1); + +for i = 1:length(spike_times) + binary_spikes(round(spike_times{i}*1000), i) = 1; +end +for i = 1:length(downsampled_stim) + start_index = (i-1) * 1000 + 1; + downsampled_stim(i) = mean(stimulus_strong()*20)) +end +