[Neuroscience] unbounded voltages in active compartmental model

dumb_founded andreajagger_8 at hotmail.com
Fri Dec 9 17:54:30 EST 2005


I am using matlab to simulate a 10 compartment model with active
Hodgkin Huxley conductances.  The problem I encounter is that the
voltages in all the compartments become infinite.  My alpha and beta
functions are correct as my separate Hodgkin Huxley model for one cell
works correctly.  The code for my compartmental model is below.  I
myself couldn't find anything egregious.


_______________________________________________________________________


clear

g_K=36;
g_Na=120;
g_Cl=.3;

V_Na=115;
V_K=-12;
V_Cl=10.613;

NCMPT = 10;


TMAX = 100;
DT = 0.025;
t = 0:DT:TMAX;


%
% specific resistances and capacitances
%
spRM = 10.0; 	% kOhm-cm^2
spCM = 1.0;	% uF/cm^2
spRA = 0.1;	% kOhm-cm

xc = 0.05 + 0.1*(0:NCMPT-1);



ISTEP = 10 ;

radius = [ 2, 2, 2, 2, 2, 2, 2, 2, 2, 2;
           2, 2, 2, 2, 2, 1, 1, 1, 1, 1;
           2, 2, 2, 2, 2, 4, 4, 4, 4, 4;
           1, 1, 1, 1, 1, 2, 2, 2, 2, 2;
           4, 4, 4, 4, 4, 2, 2, 2, 2, 2];
radius = radius*1e-4;	% convert um to cm


for ca = 1:1

    r = radius(ca,:);

    len = .1 * sqrt(spRM/(2*spRA))*sqrt(r);		% make each cmpt 0.1
lambda
    sa = 2*pi*r.*len;					% surface area
    xa = pi*r.^2;					% cross-sectional area

    RM = spRM./sa;					% compute absolute RM, CM, RA
    CM = spCM.*sa;
    RA = (spRA*len)./xa;


    v(1,1:NCMPT) = 0;					% initialization

    n(1,1:NCMPT)=alpha_n(v(1,1))/(alpha_n(v(1,1))+beta_n(v(1,1)));
    m(1,1:NCMPT)=alpha_m(v(1,1))/(alpha_m(v(1,1))+beta_m(v(1,1)));
    h(1,1:NCMPT)=alpha_h(v(1,1))/(alpha_h(v(1,1))+beta_h(v(1,1)));

    for i=2:length(t)                               	% MAIN LOOP

	 % first compartment
	 vdot(1) = ISTEP - v(i-1,1)/RM(1)...
                    - (v(i-1,1)-v(i-1,2))/((RA(1)+RA(2))/2)...

-g_K*n(i-1,1)^4*(v(i-1,1)-V_K)-g_Na*m(i-1,1)^3*h(i-1,1)*(v(i-1,1)-V_Na)-g_Cl*(v(i-1,1)-V_Cl);


	 % middle compartments
	 for c=2:NCMPT-1;
	 vdot(c) = -v(i-1,c)/RM(c) ...
		   -(v(i-1,c)-v(i-1,c-1))./((RA(c)+RA(c-1))/2) ...
                   -(v(i-1,c)-v(i-1,c+1))./((RA(c)+RA(c+1))/2)...

-g_K*n(i-1,c)^4*(v(i-1,c)-V_K)-g_Na*m(i-1,c)^3*h(i-1,c)*(v(i-1,c)-V_Na)-g_Cl*(v(i-1,c)-V_Cl);

     end


	 % last compartment
	 vdot(NCMPT) = -v(i-1,NCMPT)/RM(NCMPT) ...

-(v(i-1,NCMPT)-v(i-1,NCMPT-1))/((RA(NCMPT)+RA(NCMPT-1))/2)...

-g_K*n(i-1,NCMPT)^4*(v(i-1,NCMPT)-V_K)-g_Na*m(i-1,NCMPT)^3*h(i-1,NCMPT)*(v(i-1,NCMPT)-V_Na)-g_Cl*(v(i-1,NCMPT)-V_Cl);



	 v(i,:) = v(i-1,:) + (vdot./CM)*DT;		% Euler integration


     n_deriv=alpha_n(v(i,1))*(1-n(i-1,1))-beta_n(v(i,1))*n(i-1,1);
     m_deriv=alpha_m(v(i,1))*(1-m(i-1,1))-beta_m(v(i,1))*m(i-1,1);
     h_deriv=alpha_h(v(i,1))*(1-h(i-1,1))-beta_n(v(i,1))*h(i-1,1);
     n(i,1)=n(i-1,1)+n_deriv*DT;
     m(i,1)=m(i-1,1)+m_deriv*DT;
     h(i,1)=h(i-1,1)+h_deriv*DT;

     for c=2:NCMPT-1
     n_deriv=alpha_n(v(i,c))*(1-n(i-1,c))-beta_n(v(i,c))*n(i-1,c);
     m_deriv=alpha_m(v(i,c))*(1-m(i-1,c))-beta_m(v(i,c))*m(i-1,c);
     h_deriv=alpha_h(v(i,c))*(1-h(i-1,c))-beta_n(v(i,c))*h(i-1,c);
     n(i,c)=n(i-1,c)+n_deriv*DT;
     m(i,c)=m(i-1,c)+m_deriv*DT;
     h(i,c)=h(i-1,c)+h_deriv*DT;
     end


n_deriv=alpha_n(v(i,NCMPT))*(1-n(i-1,NCMPT))-beta_n(v(i,NCMPT))*n(i-1,NCMPT);

m_deriv=alpha_m(v(i,NCMPT))*(1-m(i-1,NCMPT))-beta_m(v(i,NCMPT))*m(i-1,NCMPT);

h_deriv=alpha_h(v(i,NCMPT))*(1-h(i-1,NCMPT))-beta_n(v(i,NCMPT))*h(i-1,NCMPT);
     n(i,NCMPT)=n(i-1,NCMPT)+n_deriv*DT;
     m(i,NCMPT)=m(i-1,NCMPT)+m_deriv*DT;
     h(i,NCMPT)=h(i-1,NCMPT)+h_deriv*DT;

     
     
    end

   
    plot(t, v);
   
end



More information about the Neur-sci mailing list