Need help in plotting current and conductance

2 vues (au cours des 30 derniers jours)
Haya Ali
Haya Ali le 3 Août 2022
Réponse apportée : KSSV le 3 Août 2022
Please help me to plot gCa, gK,ICa,IK. Below is my code
%% Forward Euler Method
close all; clear all; clc;
%Constants set for all Methods
Cm=20; % Membrane Capcitance uF/cm^2
dt=0.01; % Time Step ms
t=0:dt:200; %Time Array ms
I=90; %External Current Applied
ECa=120; % mv Na reversal potential
EK=-84; % mv K reversal potential
El=-60; % mv Leakage reversal potential
gbarCa=4.4; % mS/cm^2 Na conductance
gbarK=8; % mS/cm^2 K conductance
gbarl=2; % mS/cm^2 Leakage conductance
V(1)=-65; % Initial Membrane voltage
V1 = -1.2;
V2 = 18;
w(1)=aw(V(1))/(aw(V(1))+bw(V(1))); % Initial w-value
for i=1:length(t)-1
%Euler method to find the next m/n/h value
minf(i+1) = 1/2*(1 + tanh((V(i)-V1)/V2));
w(i+1)=w(i)+dt*((aw(V(i))*(1-w(i)))-(bw(V(i))*w(i)));
gCa=gbarCa*minf(i);
gK=gbarK*w(i);
gl=gbarl;
ICa=gCa*(V(i)-ECa);
IK=gK*(V(i)-EK);
Il=gl*(V(i)-El);
%Euler method to find the next voltage value
V(i+1)=V(i)+(dt)*((1/Cm)*(I-(ICa+IK+Il)));
end
%Store variables for graphing later
FE=V;
FEw=w;
%Plot the functions
figure
plot(t,FE);
legend('Forward Euler');
xlabel('Time (ms)');
ylabel('Voltage (mV)');
title('Voltage Change for Hodgkin-Huxley Model');
figure
plot(t,FEw,'b',t,minf,'g');
ylabel('Gaiting Variables')
xlabel('Time (ms)')
legend('w Euler','minf Euler');
figure('Name', 'Conductance')
plot(gCa, 'b', t, gK, 'g', 'LineWidth', 2)
legend('Action Potential', 'Ca+ Conductance', 'K+ Conductance')
xlabel('Time (ms)')
ylabel('Voltage (mV)')
title('Conduction of K+ and Ca+')
figure('Name', 'Currents')
plot(t,ICa, 'r',t,IK, 'b', 'LineWidth', 2)
legend('ICa+', 'IK+')
xlabel('Time (ms)')
ylabel('Current')
title ('Currents')
%%HH Function
function dydt = HH(V,y)
% Constants
ECa=120; % mv Na reversal potential
EK=-84; % mv K reversal potential
El=-60; % mv Leakage reversal potential
gbarCa=4.4; % mS/cm^2 Na conductance
gbarK=8; % mS/cm^2 K conductance
gbarl=2 % mS/cm^2 Leakage conductance
I=0; %External Current Applied
Cm=20; % Membrane Capcitance uF/cm^2
%V1 = -1.2;
%V2 = 18;
%minf = 1/2*(1 + tanh((V-V1)/V2));
% Values set to equal input values
V = y(1);
w = y(2);
gCa=gbarCa*minf;
gK=gbarK*w;
gl=gbarl;
ICa=gCa*(V-ECa);
IK=gK*(V-EK);
Il=gl*(V-El);
%Hodgkin-Huxley Model Equation
dydt = [((1/Cm)*(I-(ICa+IK+Il))); aw(V)*(1-w)-bw(V)*w];
end
%%AlphaBeta Function
function a=aw(v) %Alpha for Variable w
V3 = 2; % mV
V4 = 30; % mV
phi = 0.04; %1/ms %Rate scaling parameter
a=1/2*phi* cosh((v-V3)/(2*V4))*(1 + tanh((v-V3)/V4));
end
function b=bw(v) %Beta for variable w
V3 = 2; % mV
V4 = 30; % mV
phi = 0.04; %1/ms %Rate scaling parameter
b=1/2*phi* cosh((v-V3)/(2*V4))*(1 - tanh((v-V3)/V4));
end

Réponse acceptée

KSSV
KSSV le 3 Août 2022
You need to initialize the arrays before loop for speed and to store the values into an array.
NOTE: Check the last values of IC,IK, gCa and gK. As of now they are zeros.
clc; clear all ;
%% Forward Euler Method
close all; clear all; clc;
%Constants set for all Methods
Cm=20; % Membrane Capcitance uF/cm^2
dt=0.01; % Time Step ms
t=0:dt:200; %Time Array ms
I=90; %External Current Applied
ECa=120; % mv Na reversal potential
EK=-84; % mv K reversal potential
El=-60; % mv Leakage reversal potential
gbarCa=4.4; % mS/cm^2 Na conductance
gbarK=8; % mS/cm^2 K conductance
gbarl=2; % mS/cm^2 Leakage conductance
% Initialize
V = zeros(size(t)) ;
minf = zeros(size(t)) ;
w = zeros(size(t)) ;
gCa = zeros(size(t)) ;
gK = zeros(size(t)) ;
ICa = zeros(size(t)) ;
IK = zeros(size(t)) ;
V(1)=-65; % Initial Membrane voltage
V1 = -1.2;
V2 = 18;
w(1)=aw(V(1))/(aw(V(1))+bw(V(1))); % Initial w-value
for i=1:length(t)-1
%Euler method to find the next m/n/h value
minf(i+1) = 1/2*(1 + tanh((V(i)-V1)/V2));
w(i+1)=w(i)+dt*((aw(V(i))*(1-w(i)))-(bw(V(i))*w(i)));
gCa(i)=gbarCa*minf(i);
gK(i)=gbarK*w(i);
gl=gbarl;
ICa(i)=gCa(i)*(V(i)-ECa);
IK(i)=gK(i)*(V(i)-EK);
Il=gl*(V(i)-El);
%Euler method to find the next voltage value
V(i+1)=V(i)+(dt)*((1/Cm)*(I-(ICa(i)+IK(i)+Il)));
end
%Store variables for graphing later
FE=V;
FEw=w;
%Plot the functions
figure
plot(t,FE);
legend('Forward Euler');
xlabel('Time (ms)');
ylabel('Voltage (mV)');
title('Voltage Change for Hodgkin-Huxley Model');
figure
plot(t,FEw,'b',t,minf,'g');
ylabel('Gaiting Variables')
xlabel('Time (ms)')
legend('w Euler','minf Euler');
figure('Name', 'Conductance')
plot(t,gCa, 'b', t, gK, 'g', 'LineWidth', 2)
legend('Action Potential', 'Ca+ Conductance', 'K+ Conductance')
Warning: Ignoring extra legend entries.
xlabel('Time (ms)')
ylabel('Voltage (mV)')
title('Conduction of K+ and Ca+')
figure('Name', 'Currents')
plot(t,ICa, 'r',t,IK, 'b', 'LineWidth', 2)
legend('ICa+', 'IK+')
xlabel('Time (ms)')
ylabel('Current')
title ('Currents')
%%HH Function
%%AlphaBeta Function
function a=aw(v) %Alpha for Variable w
V3 = 2; % mV
V4 = 30; % mV
phi = 0.04; %1/ms %Rate scaling parameter
a=1/2*phi* cosh((v-V3)/(2*V4))*(1 + tanh((v-V3)/V4));
end
function b=bw(v) %Beta for variable w
V3 = 2; % mV
V4 = 30; % mV
phi = 0.04; %1/ms %Rate scaling parameter
b=1/2*phi* cosh((v-V3)/(2*V4))*(1 - tanh((v-V3)/V4));
end

Plus de réponses (0)

Catégories

En savoir plus sur Symbolic Math Toolbox dans Help Center et File Exchange

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by