动作电位的模拟程序(根据HH模型)

%------------------------------------------------------------------------%%Hodgekin-Huxleymodelofmembraneactionpotential%ref:Aquantitativedescriptionofmembranecurrentanditsapplication%toconductionandexcitationinnerve(1952,J.Physiol.,117,500-544)%Theequationswererewrittensothatoutsideofmembraneisreferencedas%zeropotential.Nopropagation.%PracticedbyWangXin-yu,inNanKaiUniversity,Dec20,2010%------------------------------------------------------------------------%clearall;GmaxNa=120;%mS/cm^2,maxconductanceEkNa=50;%mV,NerstpotentialGmaxK=36;%mS/cm^2,maxconductanceEkK=-80;%mV,NerstpotentialGvCl=0.3;%mS/cm^2,constantconductanceEkCl=-56;%mV,NerstpotentialC=1;%uF/cm^2,capacitanceofmembrane%----------initial-------------%V=-67.3;%initialmembranepotential%m=0.08;%@-48mV%n=0.32;%@-48mV%h=0.59;%@-48mV%m=0.065;%@-57.6mV%n=0.32;%@-57.6mV%h=0.59;%@-57.6mVEr=-67.3;%mV,restpotentialVr=V-Er;%potentialrelativetoEralpha_n=0.01*(-Vr+10)/(exp((-Vr+10)/10)-1);beta_n=0.125*exp(Vr/80);alpha_m=0.1*(-Vr+25)/(exp((-Vr+25)/10)-1);beta_m=4*exp(-Vr/18);alpha_h=0.07*exp(-Vr/20);beta_h=1/(exp((-Vr+30)/10)+1);m=alpha_m/(alpha_m+beta_m);n=alpha_n/(alpha_n+beta_n);h=alpha_h/(alpha_h+beta_h);GvNa=GmaxNa*m^3*h;GvK=GmaxK*n^4;I_Na=GvNa*(V-EkNa);%uAI_K=GvK*(V-EkK);%uAI_Cl=GvCl*(V-EkCl);%uA%Er=(GvNa*EkNa+GvK*EkK+GvCl*EkCl)/(GvNa+GvK+GvCl);dt=0.001;%msec,timestep%----------initialend--------------%num_idx=0;%-----------------------------------------beginfor-----------------------------------------%fortime=0:dt:4Vr=V-Er;alpha_n=0.01*(-Vr+10)/(exp((-Vr+10)/10)-1);beta_n=0.125*exp(-Vr/80);alpha_m=0.1*(-Vr+25)/(exp((-Vr+25)/10)-1);beta_m=4*exp(-Vr/18);alpha_h=0.07*exp(-Vr/20);beta_h=1/(exp((-Vr+30)/10)+1);GvNa=GmaxNa*m^3*h;GvK=GmaxK*n^4;I_Na=GvNa*(V-EkNa);%uAI_K=GvK*(V-EkK);%uAI_Cl=GvCl*(V-EkCl);%-----------------------------------changesafterdt---------------------------------%dm=(alpha_m*(1-m)-beta_m*m)*dt;dn=(alpha_n*(1-n)-beta_n*n)*dt;dh=(alpha_h*(1-h)-beta_h*h)*dt;m=m+dm;n=n+dn;h=h+dh;%-------------------injectedcurrent---------------------%iftime0.1Icm=100;%uA/cm^2elseIcm=0;end%-------------------injectedcurrent---------------------%dV=(Icm-(I_Na+I_K+I_Cl))*dt/C;%mVV=V+dV;%-----------------------------------changesafterdt---------------------------------%num_idx=num_idx+1;disp_V(num_idx)=V;disp_time(num_idx)=time;end%-----------------------------------------endfor--------------------------------------------%plot(disp_time,disp_V);是这样的：“%----------initial-------------%”之前的那几个参数如果设置的不好，形状就不对；初始的Er是可以通过自洽算出来的，图省事就直接给了，但初始V和m,n,h的配合是根据你手里那本书读图读出来的，差异不大。但如果事先不知道（没有别人的结果的话），这几个参数也是问题。所以，问题就是这些个参数如何在不知道的情况下猜出来。