added solution for fano time project

This commit is contained in:
Jan Benda 2017-01-20 18:55:35 +01:00
parent bd055b12e5
commit 52b7d39712
7 changed files with 202 additions and 0 deletions

View File

@ -0,0 +1,11 @@
function [counts, cbins] = counthist(spikes, tmin, tmax, T, cmax)
tbins = tmin+T/2:T:tmax;
cbins = 0.5:cmax;
counts = zeros(1, length(cbins));
for k = 1:length(spikes)
times = spikes{k};
n = hist(times((times>=tmin)&(times<=tmax)), tbins);
counts = counts + hist(n, cbins);
end
counts = counts / sum(counts);
end

View File

@ -0,0 +1,23 @@
function [d, thresholds, true1s, false1s, true2s, false2s, pratio] = discriminability(spikes1, spikes2, tmax, T, cmax)
[c1, b1] = counthist(spikes1, 0.0, tmax, T, cmax);
[c2, b2] = counthist(spikes2, 0.0, tmax, T, cmax);
thresholds = 0:cmax;
true1s = zeros(length(thresholds), 1);
true2s = zeros(length(thresholds), 1);
false1s = zeros(length(thresholds), 1);
false2s = zeros(length(thresholds), 1);
for k = 1:length(thresholds)
th = thresholds(k);
t1 = sum(c1(b1<=th));
f1 = sum(c1(b1>th));
t2 = sum(c2(b2>=th));
f2 = sum(c2(b2<th));
true1s(k) = t1;
true2s(k) = t2;
false1s(k) = f1;
false2s(k) = f2;
end
%pratio = (true1s + true2s)./(false1s+false2s);
pratio = (true1s + true2s)/2;
d = max(pratio);
end

View File

@ -0,0 +1,71 @@
%% general settings for the model neuron:
trials = 10;
tmax = 50.0;
D = 0.01;
%% generate and plot spiketrains for two inputs:
I1 = 14.0;
I2 = 15.0;
spikes1 = lifspikes(trials, I1, tmax, D);
spikes2 = lifspikes(trials, I2, tmax, D);
subplot(1, 2, 1);
tmin = 10.0;
spikeraster(spikes1, tmin, tmin+2.0);
title(sprintf('I_1=%g', I1))
subplot(1, 2, 2);
spikeraster(spikes2, tmin, tmin+2.0);
title(sprintf('I_2=%g', I2))
%savefigpdf(gcf(), 'spikeraster.pdf')
%% spike count histograms:
Ts = [0.01 0.1 0.5 1.0];
cmax = 100;
figure()
for k = 1:length(Ts)
T = Ts(k);
[c1, b1] = counthist(spikes1, 0.0, tmax, T, cmax);
[c2, b2] = counthist(spikes2, 0.0, tmax, T, cmax);
subplot(2, 2, k)
bar(b1, c1, 'r');
hold on;
bar(b2, c2, 'b');
xlim([0 cmax])
title(sprintf('T=%gms', 1000.0*T))
hold off;
end
%% discrimination measure:
T = 0.1;
cmax = 20;
[d, thresholds, true1s, false1s, true2s, false2s, pratio] = discriminability(spikes1, spikes2, tmax, T, cmax);
figure()
subplot(1, 3, 1);
plot(thresholds, true1s, 'b');
hold on;
plot(thresholds, true2s, 'b');
plot(thresholds, false1s, 'r');
plot(thresholds, false2s, 'r');
hold off;
% Ratio:
subplot(1, 3, 2);
fprintf('discriminability = %g\n', d);
plot(thresholds, pratio);
% ROC:
subplot(1, 3, 3);
plot(false2s, true1s);
%% discriminability:
Ts = 0.01:0.01:1.0;
cmax = 100;
ds = zeros(length(Ts), 1)
for k = 1:length(Ts)
T = Ts(k);
[c1, b1] = counthist(spikes1, 0.0, tmax, T, cmax);
[c2, b2] = counthist(spikes2, 0.0, tmax, T, cmax);
[d, thresholds, true1s, false1s, true2s, false2s, pratio] = discriminability(spikes1, spikes2, tmax, T, cmax);
ds(k) = d;
end
figure()
plot(Ts, ds)

View File

@ -0,0 +1,39 @@
function spikes = lifspikes(trials, input, tmaxdt, D)
% Generate spike times of a leaky integrate-and-fire neuron.
% trials: the number of trials to be generated
% input: the stimulus either as a single value or as a vector
% tmaxdt: in case of a single value stimulus the duration of a trial
% in case of a vector as a stimulus the time step
% D: the strength of additive white noise
tau = 0.01;
if nargin < 4
D = 1e0;
end
vreset = 0.0;
vthresh = 10.0;
dt = 5e-5;
if max(size(input)) == 1
input = input * ones(ceil(tmaxdt/dt), 1);
else
dt = tmaxdt;
end
spikes = cell(trials, 1);
for k=1:trials
times = [];
j = 1;
v = vreset + (vthresh-vreset)*rand(1);
noise = sqrt(2.0*D)*randn(length(input), 1)/sqrt(dt);
for i=1:length(noise)
v = v + (- v + noise(i) + input(i))*dt/tau;
if v >= vthresh
v = vreset;
spiketime = i*dt;
times(j) = spiketime;
j = j + 1;
end
end
spikes{k} = times;
end
end

View File

@ -0,0 +1,28 @@
function savefigpdf(fig, name, width, height)
% Saves figure fig in pdf file name.pdf with appropriately set page size
% and fonts
% default width:
if nargin < 3
width = 11.7;
end
% default height:
if nargin < 4
height = 9.0;
end
% paper:
set(fig, 'PaperUnits', 'centimeters');
set(fig, 'PaperSize', [width height]);
set(fig, 'PaperPosition', [0.0 0.0 width height]);
set(fig, 'Color', 'white')
% font:
set(findall(fig, 'type', 'axes'), 'FontSize', 12)
set(findall(fig, 'type', 'text'), 'FontSize', 12)
% save:
saveas(fig, name, 'pdf')
end

View File

@ -0,0 +1,30 @@
function spikeraster(spikes, tmin, tmax)
% Display a spike raster of the spike times given in spikes.
%
% spikeraster(spikes, tmax)
% spikes: a cell array of vectors of spike times in seconds
% tmin: plot spike raster starting at tmin seconds
% tmax: plot spike raster upto tmax seconds
ntrials = length(spikes);
for k = 1:ntrials
times = spikes{k};
times = times((times>=tmin) & (times<=tmax));
if tmax < 1.5
times = 1000.0*times; % conversion to ms
end
for i = 1:length( times )
line([times(i) times(i)],[k-0.4 k+0.4], 'Color', 'k');
end
end
if (tmax-tmin) < 1.5
xlabel('Time [ms]');
xlim([1000.0*tmin 1000.0*tmax]);
else
xlabel('Time [s]');
xlim([tmin tmax]);
end
ylabel('Trials');
ylim([0.3 ntrials+0.7 ]);
end