Solving FitzHugh-Nagumo equations using ODE45

Write a program to solve the FitzHugh-Nagumo equations for a single cell (i.e., without spatial coupling).

du/dt = c1u ( u − a)(1 − u) − c2uv +stim

dv/ dt = b (u − v)

where

a=0.13

b=0.013

c1=0.26

c2=0.1

stim is a stimulus current that can be applied for a short time at the beginning of the simulation.

u represents membrane potential and ranges from 0 (rest) to 1 (excited). v is a recovery variable in the same range. t is time in milliseconds.

How do you use MATLAB's ode45() function to integrate the system of differential equations? Input to the program should be the duration of the simulation; initial values for u, v, and t; the strength of the stimulus, and the time for which it is applied (typically a few ms). It;s output should include vectors for t, u and v.

[darova]

Once you mentioned ode45. Did you read help?

[Kate Heinzman]

Yes but I'm still confused as to how to write the overall function with the specified input and output arguments stated above.

Answer 1

[darova]

Here is my achievement. I don't get about stim. Can you explain more?

%             du/dt = c1u ( u  − a)(1 − u) − c2uv +stim
%             dv/ dt  = b (u  − v)  
a = 0.13;
b = 0.013;
c1 = 0.26;
c2 = 0.1;
% let y(1) = u; y(2) = v
F = @(t,y)  [c1*y(1)*(y(1)-a)*(1-y(1)-c2*y(1)*y(2)+stim)
            b*(y(1)-y(2))];
tspan = [0 2];              % time
y0 = [1 2];                 % u0 = 1; v0 = 2;
[t,y] = ode45(F,tspan,y0);
plot(t,y)
legend('u(t)','v(t)')

[Kate Heinzman]

Incorporating stim is what I am confused about. Here is the code I have so far, but it's giving me errors

function [t,u,v] = fhn0D(tf,u0,v0,t0,stim,tstim)

% Solve the FitzHugh-Nagumo equations for a single cell (i.e., without

% spatial coupling)

% du/dt = (c1*u)*(u − a)*(1 − u) − (c2*u*v) + stim

% dv/d t = b(u − v)

% stim is a stimulus current that can be applied for a short time at

% the beginning of the simulation

% u represents membrane potential and ranges from 0 (rest) to 1 (excited)

% v is a recovery variable in the same range

% t is time in milliseconds

% INPUT: tf duration of simulation

% u0 initial value of u

% v0 initial value of v

% t0 initial value of t

% stim strength of stimulus

% tstim time of applied stimulus

%

% OUTPUT: t time vector (ms)

% u membrane potential vector

% v recovery vector

[t,U] = ode45(@fhn,[t0 tf],[u0; v0]);

% First and second column of U correspond to u and v respectively

u = U(:,1);

v = U(:,2);

% Plot u and v vs t

plot(t,u,t,v)

title('Solution of FitzHugh-Nagumo Equations for a Single Cell with ODE45');

xlabel('Time t');

ylabel('Solution U');

legend('u','v');

end

function dUdt = fhn(t,U)

% FitzHugh-Nagumo equations defined to be plugged into ode45

%

% INPUT: t time

% U vector that holds u and v

%

% OUTPUT: dUdt vector containing solutions to derivatives of u and v

% U(1) is u and U(2) is v

% Other needed variables

a = 0.13;

b = 0.013;

c1 = 0.26;

c2 = 0.1;

% Preallocate the output vector

dUdt = zeros(2,1);

% Create a two element vector that holds the derivative equations of u and

% v

dUdt = [(c1*U(1))*(U(1)-a)*(1-U(1))-(c2*U(1)*U(2))+ stim; b*(U(1)-U(2))];

end

[darova]

• stim is a stimulus current that can be applied for a short time at the beginning of the simulation.

As i understand correctly:

stim = 2;
tstim = 0.5;
F = @(t,y)  [c1*y(1)*(y(1)-a)*(1-y(1)-c2*y(1)*y(2)+stim*(t<tstim))
            b*(y(1)-y(2))];

the same as

if t < tstim
    stim = 2;
else
    stim = 0;
end

[Kate Heinzman]

I understand what you're saying but since stim is an input argument in the first function in my code it gives me an error when I try to use it in my second function. How do I pass stim from the first function to the second function?

[darova]

Just add it to input arguments?

function dUdt = fhn(t,U,stim)

Then call function

[t,U] = ode45(@(t,U)fhn(t,U,stim),[t0 tf],[u0; v0]);