MATLAB Answers

Finding membrane time constant from a graph

1 view (last 30 days)
I wrote a script to numerically calculate membrane potential produced by intracellular current injection in a cell model using euler integration method to compare it to the analytical membrane potential. I'm not sure how to calculate the membrane time constant from the numerical Vm(t)as a time to reach ~ 0.632 of the maximal value.
% Model parameters
a = 0.001; % cell radius in cm
Rm = 5000; % specific membrane resistance in ohm*cm^2
Cm = 0.001; % specific membrane capacitance in mF/cm^2
I0 = 0.1; % current strength in mA
T = 50; % integration duration in ms
t = 0.02; % integration step in ms
Tm = Rm * Cm; % time constant in ms
num_iter = T / t; % number of interations of time loop
% Preallocate Vm_num (numerical Vm), Vm_ana (analytical Vm), and time
% vectors for the for loop
Vm_num = zeros(1,num_iter);
Vm_ana = zeros(1,num_iter);
time = zeros(1,num_iter);
% When calculating Rm * I0 in the Euler method, Rm is in terms of the area
% of the cell so need to multiply it by the surface area to get rid of the
% cm^2 in the units. So Rm * I0 is in terms of mV instead of mV * cm^2
% Use surface area bc that is where the shock is being applied
% Surface area of a circle is 4*pi*r^2
RM = Rm / (4*pi*a^2);
for j = 1:num_iter
% initial condition: Vm(1) = 0 so start at 2nd index
Vm_num(j+1) = Vm_num(j) + (t/Tm)*((RM*I0)-Vm_num(j));
Vm_ana(j+1) = (RM*I0)*(1 - exp(-time(j)/Tm));
time(j+1) = time(j) + t;
end
% Plot numerical and analytical Vm on the same plot
plot(time,Vm_num);
hold on
plot(time,Vm_ana);
title('Comparing Analytical and Numerical Membrane Potentials over Time');
xlabel('Time(ms)');
ylabel('Membrane Potential (mV)');
legend('Numerical Vm','Analytical Vm');

  0 Comments

Sign in to comment.

Accepted Answer

Star Strider
Star Strider on 5 Apr 2020
Try adding this line to your code (after the loop):
tau_Vm_num = interp1(Vm_num, time, 0.632*max(Vm_num));
so the entire code is now:
% Model parameters
a = 0.001; % cell radius in cm
Rm = 5000; % specific membrane resistance in ohm*cm^2
Cm = 0.001; % specific membrane capacitance in mF/cm^2
I0 = 0.1; % current strength in mA
T = 50; % integration duration in ms
t = 0.02; % integration step in ms
Tm = Rm * Cm; % time constant in ms
num_iter = T / t; % number of interations of time loop
% Preallocate Vm_num (numerical Vm), Vm_ana (analytical Vm), and time
% vectors for the for loop
Vm_num = zeros(1,num_iter);
Vm_ana = zeros(1,num_iter);
time = zeros(1,num_iter);
% When calculating Rm * I0 in the Euler method, Rm is in terms of the area
% of the cell so need to multiply it by the surface area to get rid of the
% cm^2 in the units. So Rm * I0 is in terms of mV instead of mV * cm^2
% Use surface area bc that is where the shock is being applied
% Surface area of a circle is 4*pi*r^2
RM = Rm / (4*pi*a^2);
for j = 1:num_iter
% initial condition: Vm(1) = 0 so start at 2nd index
Vm_num(j+1) = Vm_num(j) + (t/Tm)*((RM*I0)-Vm_num(j));
Vm_ana(j+1) = (RM*I0)*(1 - exp(-time(j)/Tm));
time(j+1) = time(j) + t;
end
tau_Vm_num = interp1(Vm_num, time, 0.632*max(Vm_num));
% Plot numerical and analytical Vm on the same plot
plot(time,Vm_num);
hold on
plot(time,Vm_ana);
plot(tau_Vm_num, 0.632*max(Vm_num), 'sr')
title('Comparing Analytical and Numerical Membrane Potentials over Time');
xlabel('Time(ms)');
ylabel('Membrane Potential (mV)');
legend('Numerical Vm','Analytical Vm','\tau');
with ‘tau_Vm_num’ being the time constant, and the square marker in the plot denoting the time constant as you have defined it.

More Answers (0)

Sign in to answer this question.