EE746 Neuromorphic Engineering Assignment 1 Solved.

       Assignment 1                                                                                         Neuromorphic Engineering
Note: In all the problems, be mindful of the units of various quantities and the sign conventions for currents and voltages.

Your submission must include:

1. A .pdf or .doc file clearly documenting your code, figures, and results.

2. Your MATLAB code. The code must be saved in plain text files that can be immediately runin MATLAB; include a file called main.m that runs the functions you implemented and generates the figures described in your assignment write-up. You will lose credit if your code is absent or cannot be run.

Store your write-up and code in a single directory named hw1 yourID (for example hw1 12D423222) and submit it in moodle.

Do not wait till the last minute to start the Assignment as you will require at least at least 8 hours to complete all the questions.

Late submission policy:

Before the solution key is uploaded in moodle: If your original score is S and you submitted the HW X hours after the deadline, your score will be S exp(−X/24).

After the solution key is uploaded in moodle: 0 credit.

Spiking Neuron Models
In this assignment, we will learn how to model the activity of spiking neurons starting with the simplest model. We will also implement the Hodgkin-Huxley neuron model and use it to determine an estimate of the energy cost associated with a spike.

Problem 1: Leaky Integrate and Fire Model
The dynamics of the membrane potential V (t), in the LIF neuron model is given by the equation

                                                      )                                                (1)

C is the membrane capacitance, gL is the leak conductance and EL is the leak reversal potential. Iapp is the externally applied current and is assumed to be positive for current flowing into the cell.

When V ≥ VT , a spike is issued and V is reset to EL (We will write this as V (t) −→ EL). (Assume that C= 300 pF, gL = 30 nS, VT = 20 mV and EL = −70 mV).

(a)  Assume that a constant current I0 is applied to the neuron. Write an expression for the steady state value of the membrane potential. Hence, determine the minimum value of the steady state current, Ic necessary to initiate a spike. 2 points

(b)  In order to simulate the behavior of a set of N neurons, it is useful to define a N × 1 column vector to the store the membrane potentials of the N neurons. Write a program to solve the equivalent difference equation numerically using Runge-Kutta second order method for a set of N neurons driven by external current input. (You should not use a FOR loop to calculate the potential of the N neurons, rather, the potential of all the neurons should be calculated in one step.)

Assume that the input to the program is a N ×M column vector, representing the input current for the N neurons for M time-intervals where M = T/∆t. The output of your program should be a N × M matrix storing the values of the membrane potential for the N neurons, for the M time-intervals.       8 points

(c)  We would like to now use this framework to study the dynamics of LIF neurons. Assumethat you have a population of 10 identical neurons, with each neuron receiving a constant current. Let the magnitude of the input current for the kth neuron be given by the expression

                                                                      Iapp,k = (1 + kα)Ic                                                                                                   (2)

where α = 0.1. (In this example, the current does not vary with time, so, all the values across any row of the input current matrix is a constant). Plot the membrane potential for neurons 2,4,6 and 8 from t = 0 to 500 ms. (Assume ∆t = 0.1 ms and at t = 0, the neuron is in steady-state with

Iapp(t) = 0).
6 points
(d) Plot the average time interval between spikes from (c) as a function of Iapp,k.
4 points
Problem 2: Izhikevich Model
The dynamics of the membrane potential V (t), in the Izhikevich neuron model is given by the equations

                                                                                                                                   )                                (3)

                                                                                                )]                                                                (4)

When V (t) ≥ vpeak, V (t) −→ c and U(t) −→ U(t) + d.

By varying the parameters C,Er,Et,kz,a,b,c and d, a variety of neuronal behaviors can be modeled.

 
C(pF)
kz(µS/V) Er(mV)
Et(mV)
a (KHz) b(nS)
c(mV)
d(pA)
vpeak(mV)
RS
100
0.7
−60
−40
0.03
−2
−50
100
35
IB
150
1.2
−75
−45
0.01
+5
−56
130
50
CH
50
1.5
−60
−40
0.03
+1
−40
150
25
(a) What are the steady state values of V and U for Iapp = 0?
2 points
 
(b) Write the equivalent difference equations for (3) and (4).
2 points
 
(c) Write a program to solve the equivalent difference equation for a set of N Izhikevich model neurons using Runge-Kutta fourth order method. The neuron type should be a parameter in your function call, for each of the N neurons. You may use ∆t = 0.1 ms and plot the response of the three neurons above from t = 0 to 500 ms, for Iapp = 400,500,600 pA. 16 points

Note: Try to re-run 2(c) with larger values of ∆t. You will notice that the overall dynamics can change drasticially, especially for neuron CH. This is because of the inaccuracy in determining the exact time when V (t) exceeds vpeak. For a good description of this issue and how to get around it, see: Hybrid spiking models, Phil. Trans. R. Soc. A November 13, 2010 368 (1930) 5061-5070.

We will ignore these issues here, and for the sake of displaying, you may chose to artificially set the value of membrane potential just before reset to vpeak in your code.

Problem 3: Adaptive Exponential Integrate-and-Fire Model
The dynamics of the membrane potential V (t), in the AEF neuron model is given by the equations

                                                                                                                                                  )                 (5)

                                                                               )                                                                                   (6)

When V (t) ≥ 0, V (t) −→ Vr and U(t) −→ U(t) + b.

As before, by varying the parameters a variety of neuronal behaviors can be modeled.

 
C(pF)
gL(ns)
EL(mV)
VT (mV)
∆T (mV) a(nS)
τw(ms)
b(pA)
Vr(mV)
RS
200
10
−70
−50
2
2
30
0
−58
IB
130
18
−58
−50
2
4
150
120
−50
CH
200
10
−58
−50
2
2
120
100
−46
(a)  Write the equivalent difference equations for (5) and (6).           2 points

(b)  What are the steady state values of V and U for Iapp = 0? (Determine the answer numerically, such that the value of V is accurate within ±1µV). 6 points

(c)  Write a program to solve the equivalent difference equation for a set of N AEF model neurons using Euler method. The neuron type should be a parameter in your function call, for each of the N neurons. You may use ∆t = 0.1 ms and plot the response of the three neurons above from t = 0 to 500 ms, for Iapp = 250,350,450 pA.      12 points Problem 4: Spike energy based on Hodgkin-Huxley neuron model

The dynamics of the membrane potential V (t), in the Hodgkin-Huxley neuron model is given by the equations

                                                 )                                           (7)

where
 
                                                         iNa(t)    =     gNam3h(V (t) − ENa)
(8)
                                                          iK(t)    =     gKn4(V (t) − EK)
(9)
                                                            il(t)    =    gl(V (t) − El)

The variables n,m and h lie in the interval [0,1] and obey the equation
(10)
                                                                                                                    (11)

where

                                                                                              βn(t) = 0.125exp(−(V (t) + 65)/80)         (12)

                                                                                            βm(t) = 4exp(−0.0556(V (t) + 65))           (13)

                                                                                                         (14)

V (t) is assumed to be measured in mV in the above expressions.

Assume that C=1µF/cm2, ENa = 50 mV, EK = −77 mV and El = −55 mV and gNa = 120 mS/cm2, gK = 36 mS/cm2 and gl = 0.3 mS/cm2.

(a)  Solve the equivalent difference equations for (7) and (11) to determine the ion currents andthe membrane potential for a step current waveform described as follows:

Iext(t) = I0[H(t−2T)−H(t−3T)] where H(x) is the Heaviside step function defined as H(x) = 1 if x ≥ 0 and H(x) = 0 otherwise. You should set your initial conditions such that there are no spikes when there is no external input current to the neuron. Assume I0 = 15µA/cm2, ∆t = 0.01 ms, T = 30ms.     12 points

(b)  The instantaneous power dissipated (per unit area) in the various ion channels can beapproximated by the expression

                                                                  Px(t) = ix(t)(V (t) − Ex)                                                          (15)

Similarly, the power spent in charging/discharging the membrane capacitance (per unit area) can be approximated as

                                        )                                (16)

Plot the relative magnitudes of the instantaneous power dissipated in the three ion channels and in the membrane capacitor as a function of time for one cycle of the action potential. 4 points

(c)  Numerically integrate the power in (b) to determine the total energy dissipated in one cycleof the action potential for a patch of the cell membrane with area of 1 µm2 4 points.