MATH 3104/7134: Assignment 1
Prof Geoffrey Goodhill, Semester 1, 2021
This assignment is marked out of 13 but is worth 8.75% of the total mark for the course overall (6% for
MATH7134). It is due to be submitted by 4pm on 9/3/20.
What you hand in must be your own unique work. You should hand in the following:
• Fully-worked answers to the analysis questions below.
• Your Matlab code (EXTENSIVELY COMMENTED - if it’s not clear to us what the program’s doing we
can’t give you the marks).
• The results of running your program as asked for below.
A few points to remember for the Matlab component:
• Don’t hard code variables throughout your code. Set up all your variables at the start of your m-file.
• Put a title and axis labels on your plots.
• Convert all variables to SI or SI derived units.
The membrane equation
The following questions are based around the membrane equation:
τ
dVm
dt
= −Vm + Vrest +RmIinj . (1)
For explanation of terminology etc. see the lecture notes.
Analysis questions
These questions fill in some of the derivations that were skipped in the lectures. They are worth one mark
each.
1. For Vm = Vreset at t = 0 and Iinj a constant, solve the membrane equation to obtain the following
solution:
Vm = Vrest +RIinj + (Vreset − Vrest −RIinj) e−t/τ .
2. Show that the time tisi when Vm reaches Vthreshold is
tisi = τ ln
(
RIinj + Vrest − Vreset
RIinj + Vrest − Vthreshold
)
.
3. Defining a firing rate by r = 1/tisi, use ln(1+z) ≈ z for small z to show that when Iinj (again constant)
is sufficiently large,
r ≈ Vrest − Vthreshold +RIinj
τ(Vthreshold − Vreset) .
Matlab questions
4. Build an integrate-and-fire neuron from the membrane equation. Use Euler’s method1 to numerically
integrate equation 1 above. Use Vrest = −70 mV,Rm = 10 MΩ, and τ = 10 ms. Initially set Vm = Vrest.
When the membrane potential reaches Vth = −54 mV, manually cause a spike and then instantly reset
the potential to Vreset = −80 mV. Show sample voltage traces (with spikes) for a 300 ms-long current
pulse. Choose a reasonable current Iinj centred in a 500 ms-long simulation (i.e. 100 ms off, 300 ms
on, 100 ms off). [4 marks]
HINT: If your code is running slowly, make sure you’ve preallocated the size of your voltage vector
before entering the for loop. This avoids a problem arising from how Matlab handles memory. You
can do this with some code like V = zeros(n,1), where n is the number of timesteps in your for
loop. This creates a column vector of length n of zeros.
5. From your simulations determine the firing rate of the model for constant Iinj ranging from 0 to 10 nA,
and plot a graph comparing the results from your simulations with the firing rate predicted for large
Iinj calculated in question 3 above. Discuss the match. [3 marks]
6. Now include an extra current in the integrate-and-fire model to introduce spike-rate adaptation as
described in the lecture. Use τsra = 100 ms, ∆gsra = 0.3 (with rm = 1), VK = -70 mV and an initial
value of gsra of zero. Note that gsra is defined by a differential equation, so you will have to integrate it
with Euler’s method along with Vm, as well as adding ∆gsra to it at certain time points. Show a sample
voltage trace illustrating the adaptation of the interspike interval. [3 marks]
1If you have not come across this before then look it up. Note you need to make the step size small enough to ensure numerical
stability.