Leaky Integrate-and-Fire Model (LIF Model) matlab code

Set up the default plotting parameters

set(0,'DefaultLineLineWidth',2,...
    'DefaultLineMarkerSize',8, ...
    'DefaultAxesLineWidth',2, ...
    'DefaultAxesFontSize',14,...
    'DefaultAxesFontWeight','Bold');
figure(1)
clf;

Simulation Parameters

dt = 0.0001;                % time-step
t = 0:dt:0.5;               % vector of time-points
ton = 0.15;                 % time to begin applied current (onset)
toff = 0.35;                % time to end applied current (offset)
non = round(ton/dt);        % time-point index of current onset
noff = round(toff/dt);      % time-point index of current offset

Parameters for the LIF neuron

tau = 0.010;                % membrane time constant
E_L = -0.070;               % leak potential (also resting potential)
Vth = -0.050;               % threshold potential (to produce spike)
Vreset = -0.080;            % reset potential (post-spike)
Cm = 100e-12;               % total membrane capacitance
G_L = Cm/tau;               % total membrane conductance (leak conductance)

LIF Neuron Simulation

Iapp = 210e-12;                 % value of applied current steps
I = zeros(size(t));             % vector for current at each time-point
I(non:noff) = Iapp;             % add the applied current
V = E_L*ones(size(t));          % initialize the membrane potential vector
spikes = zeros(size(t));        % initialize a vector to record spikes

for i = 2:length(t);            % loop through all time points
    % next line: Forward Euler method to update membrane potential
    V(i) = V(i-1) + dt*(I(i) +G_L*(E_L-V(i-1)))/Cm;
    if V(i) > Vth;              % if potential is above threshold
        spikes(i) = 1;          % record the spike at that time-point
        V(i) = Vreset;          % reset the potential
    end;
end;                            % end the loop & go to next time-point

Plotting Current Against Time

figure(1)
plot(t,I*1e9,'k');              % Plot current (in nA) against time
ylabel('I_{app} (nA)')

% The following line adds a title, which uses the command "strcat" to
% combine the value of current in that trial with "nA" to make a
% complete title.
% The command "num2str" converts the number representing the applied
% current into a string, which is something that can be printed.
title(strcat(num2str(1e9*Iapp),'nA'))

set(gca,'XTick',[0 0.25 0.5])

axis([0 0.5 0 0.3])         % Sets x-axis then y-axis ranges

Plotting Membrane Potential Against Time

figure(2)
plot(t,1000*V,'k');                      % Plot membrane potential vs time
ylabel('V_m (mV)')          ㄴ     % Label y-axis
set(gca,'XTick',[0 0.25 0.5])
axis([0 0.5 1000*(Vreset-0.005) 1000*(Vth+0.005)])    % Set the axes

Plotting Spikes Against Time

figure(3)
plot(t,spikes,'k')          % Plots a line from 0 to 1 at each spike
xlabel('Time (sec)')        % Label the x-axis
ylabel('Spikes')        % Label the y-axis
set(gca,'XTick',[0 0.25 0.5])
set(gca,'YTick',[0 1])
axis([0 0.5 0 1])           % Set the ranges of x- and y-axes