[projects] added solution to isicorrelations

projects/project_isicorrelations/solution/firingrate.m

@@ -0,0 +1,11 @@
+function [ rate ] = firingrate(spikes, t0, t1)
+% compute the firing rate from spikes between t0 and t1
+
+rates = zeros(length(spikes), 1);
+for k = 1:length(spikes)
+    spiketimes = spikes{k};
+    rates(k) = length(spiketimes((spiketimes>=t0)&(spiketimes<=t1)))/(t1-t0);
+end
+rate = mean(rates);
+end
+

projects/project_isicorrelations/solution/isicorrelations.m

@@ -0,0 +1,36 @@
+trials = 5;
+tmax = 10.0;
+Dnoise = 1e-2; % noise strength
+adapttau = 0.1; % adaptation time constant in seconds
+adaptincr = 0.5; % adaptation strength
+t0=2.0;
+t1=tmax;
+maxlag = 5;
+
+taus = [0.01, 0.1, 1.0];
+colors = ['r', 'b', 'g'];
+for j = 1:length(taus)
+    adapttau = taus(j);
+    % f-I curves:
+    Is = [0:10:80];
+    rate = zeros(length(Is),1);
+    corr = zeros(length(Is),maxlag);
+    for k = 1:length(Is)
+        input = Is(k);
+        spikes = lifadaptspikes(trials, input, tmax, Dnoise, adapttau, adaptincr);
+        rate(k) = firingrate(spikes, t0, t1);
+        corr(k,:) = serialcorr(spikes, t0, t1, maxlag);
+    end
+
+    figure(1);
+    hold on
+    plot(Is, rate, colors(j));
+    hold off
+
+    figure(2);
+    hold on
+    plot(Is, corr(:,2), colors(j));
+    hold off
+end
+
+pause

projects/project_isicorrelations/solution/lifadaptspikes.m

@@ -0,0 +1,50 @@
+function spikes = lifadaptspikes( trials, input, tmaxdt, D, tauadapt, adaptincr )
+% Generate spike times of a leaky integrate-and-fire neuron
+% with an adaptation current
+% 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
+% tauadapt: adaptation time constant
+% adaptincr: adaptation strength
+    
+    tau = 0.01;
+    if nargin < 4
+        D = 1e0;
+    end
+    if nargin < 5
+        tauadapt = 0.1;
+    end
+    if nargin < 6
+        adaptincr = 1.0;
+    end
+    vreset = 0.0;
+    vthresh = 10.0;
+    dt = 1e-4;
+    
+    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();
+        a = 0.0;
+        noise = sqrt(2.0*D)*randn( length( input ), 1 )/sqrt(dt);
+        for i=1:length( noise )
+            v = v + ( - v - a + noise(i) + input(i))*dt/tau;
+            a = a + ( - a )*dt/tauadapt;
+            if v >= vthresh
+                v = vreset;
+                a = a + adaptincr/tauadapt;
+                times(j) = i*dt;
+                j = j + 1;
+            end
+        end
+        spikes{k} = times;
+    end
+end

projects/project_isicorrelations/solution/serialcorr.m

@@ -0,0 +1,15 @@
+function [ isicorrs ] = serialcorr(spikes, t0, t1, maxlag)
+% compute the serial correlations from spikes between t0 and t1
+
+ics = zeros(length(spikes), maxlag);
+for k = 1:length(spikes)
+    spiketimes = spikes{k};
+    isis = diff(spiketimes((spiketimes>=t0)&(spiketimes<=t1)));
+    if length(isis) > 2*maxlag
+        for j=1:maxlag
+            ics(k, j) = corr(isis(j:end)', isis(1:end+1-j)');
+        end
+    end
+end
+isicorrs = mean(ics, 1);
+end