/* * Mathematical Model Of A Human Atrial Action Potential, 1998 * * Model Status * * This model has been validated by Penny Noble of Oxford University * and is known to run in COR and PCEnv. * * ValidateCellML verifies this model as valid CellML, but detects * unit inconsistencies. * * Model Structure * * The ionic mechanisms underlying many important properties of * the human atrial action potential are poorly understood. Using * specific formulations of the K+, Na+ and Ca2+ currents based * on experimental data recorded from human atrial myocytes, along * with representations of pump, exchange and background currents, * Marc Courtemanche, Rafael J. Ramirez and Stanley Nattel developed * a mathematical model of the action potential (see the figure * below). * * This mathematical model builds mostly on the classical work * of Luo and Rudy (see The Luo-Rudy Ventricular Model II (dynamic), * 1994). Courtemanche et al effectively develop a working model * of the human atrial action potential from the Luo-Rudy model * which is based on guinea pig ventricular cells. Their primary * goal was to develop a useful model of the action potential from * which they could gain insights into experimental observations * made on human atrial cells and tissues and make predictions * regarding the behaviour of these cells under previously untested * conditions * * The complete original paper reference is cited below: * * Ionic mechanisms underlying human atrial action potential properties: * insights from a mathematical model, Marc Courtemanche, Rafael * J. Ramirez and Stanley Nattel, 1998, American Journal of Physiology * , 275, H301-H321. PubMed ID: 9688927 * * schematic diagram of a human atrial myocyte * * [[Image file: courtemanche_1998.png]] * * A schematic representation of currents, pumps and exchangers * included in the model. The cell includes three intracellular * compartments: cytoplasm, sarcoplasmic reticulum (SR) release * compartment [junctional SR (JSR)], and SR uptake compartment * [network SR (NSR)]. */ import nsrunit; unit conversion on; // unit millisecond predefined unit per_millisecond=1E3 second^(-1); // unit millivolt predefined unit per_millivolt=1E3 kilogram^(-1)*meter^(-2)*second^3*ampere^1; unit per_millivolt_millisecond=1E6 kilogram^(-1)*meter^(-2)*second^2*ampere^1; // unit millimolar predefined unit millimolar_per_millisecond=1E3 meter^(-3)*second^(-1)*mole^1; unit millisecond_per_millimolar=.001 meter^3*second^1*mole^(-1); unit picoA=1E-12 ampere^1; unit picoA_per_picoF=1 kilogram^1*meter^2*second^(-4)*ampere^(-1); unit nanoS_per_picoF=1E3 second^(-1); unit picoF=1E-12 kilogram^(-1)*meter^(-2)*second^4*ampere^2; unit micrometre_3=1E-18 meter^3; unit per_micrometre_3=1E18 meter^(-3); unit mm2=1E-6 meter^2; unit joule_per_mole_kelvin=1 kilogram^1*meter^2*second^(-2)*kelvin^(-1)*mole^(-1); unit coulomb_per_millimole=1E3 second^1*ampere^1*mole^(-1); math main { realDomain time millisecond; time.min=0; extern time.max; extern time.delta; real V(time) millivolt; when(time=time.min) V=-81.18; real R joule_per_mole_kelvin; R=8.3143; real T kelvin; T=310; real F coulomb_per_millimole; F=96.4867; real Cm picoF; Cm=100; real i_st(time) picoA; real i_Na(time) picoA; real i_K1(time) picoA; real i_to(time) picoA; real i_Kur(time) picoA; real i_Kr(time) picoA; real i_Ks(time) picoA; real i_Ca_L(time) picoA; real i_CaP(time) picoA; real i_NaK(time) picoA; real i_NaCa(time) picoA; real i_B_Na(time) picoA; real i_B_Ca(time) picoA; real stim_start millisecond; stim_start=50; real stim_end millisecond; stim_end=50000; real stim_period millisecond; stim_period=1000; real stim_duration millisecond; stim_duration=2; real stim_amplitude picoA; stim_amplitude=-2000; real E_Na(time) millivolt; real g_Na nanoS_per_picoF; g_Na=7.8; real Na_i(time) millimolar; when(time=time.min) Na_i=1.117e+01; real Na_o millimolar; Na_o=140; real m(time) dimensionless; when(time=time.min) m=2.908e-3; real h(time) dimensionless; when(time=time.min) h=9.649e-1; real j(time) dimensionless; when(time=time.min) j=9.775e-1; real alpha_m(time) per_millisecond; real beta_m(time) per_millisecond; real m_inf(time) dimensionless; real tau_m(time) millisecond; real alpha_h(time) per_millisecond; real beta_h(time) per_millisecond; real h_inf(time) dimensionless; real tau_h(time) millisecond; real alpha_j(time) per_millisecond; real beta_j(time) per_millisecond; real j_inf(time) dimensionless; real tau_j(time) millisecond; real E_K(time) millivolt; real g_K1 nanoS_per_picoF; g_K1=0.09; real K_o millimolar; K_o=5.4; real K_i(time) millimolar; when(time=time.min) K_i=1.39e+02; real K_Q10 dimensionless; K_Q10=3; real g_to nanoS_per_picoF; g_to=0.1652; real oa(time) dimensionless; when(time=time.min) oa=3.043e-2; real oi(time) dimensionless; when(time=time.min) oi=9.992e-1; real alpha_oa(time) per_millisecond; real beta_oa(time) per_millisecond; real tau_oa(time) millisecond; real oa_infinity(time) dimensionless; real alpha_oi(time) per_millisecond; real beta_oi(time) per_millisecond; real tau_oi(time) millisecond; real oi_infinity(time) dimensionless; real g_Kur(time) nanoS_per_picoF; real ua(time) dimensionless; when(time=time.min) ua=4.966e-3; real ui(time) dimensionless; when(time=time.min) ui=9.986e-1; real alpha_ua(time) per_millisecond; real beta_ua(time) per_millisecond; real tau_ua(time) millisecond; real ua_infinity(time) dimensionless; real alpha_ui(time) per_millisecond; real beta_ui(time) per_millisecond; real tau_ui(time) millisecond; real ui_infinity(time) dimensionless; real g_Kr nanoS_per_picoF; g_Kr=0.029411765; real xr(time) dimensionless; when(time=time.min) xr=3.296e-5; real alpha_xr(time) per_millisecond; real beta_xr(time) per_millisecond; real tau_xr(time) millisecond; real xr_infinity(time) dimensionless; real g_Ks nanoS_per_picoF; g_Ks=0.12941176; real xs(time) dimensionless; when(time=time.min) xs=1.869e-2; real alpha_xs(time) per_millisecond; real beta_xs(time) per_millisecond; real tau_xs(time) millisecond; real xs_infinity(time) dimensionless; real g_Ca_L nanoS_per_picoF; g_Ca_L=0.12375; real Ca_i(time) millimolar; when(time=time.min) Ca_i=1.013e-4; real d(time) dimensionless; when(time=time.min) d=1.367e-4; real f(time) dimensionless; when(time=time.min) f=9.996e-1; real f_Ca(time) dimensionless; when(time=time.min) f_Ca=7.755e-1; real d_infinity(time) dimensionless; real tau_d(time) millisecond; real f_infinity(time) dimensionless; real tau_f(time) millisecond; real f_Ca_infinity(time) dimensionless; real tau_f_Ca millisecond; real Km_Na_i millimolar; Km_Na_i=10; real Km_K_o millimolar; Km_K_o=1.5; real i_NaK_max picoA_per_picoF; i_NaK_max=0.59933874; real f_NaK(time) dimensionless; real sigma dimensionless; real i_B_K(time) picoA; real g_B_Na nanoS_per_picoF; g_B_Na=0.0006744375; real g_B_Ca nanoS_per_picoF; g_B_Ca=0.001131; real g_B_K nanoS_per_picoF; g_B_K=0; real E_Ca(time) millivolt; real Ca_o millimolar; Ca_o=1.8; real I_NaCa_max picoA_per_picoF; I_NaCa_max=1600; real K_mNa millimolar; K_mNa=87.5; real K_mCa millimolar; K_mCa=1.38; real K_sat dimensionless; K_sat=0.1; real gamma dimensionless; gamma=0.35; real i_CaP_max picoA_per_picoF; i_CaP_max=0.275; real i_rel(time) millimolar_per_millisecond; real Fn(time) dimensionless; real K_rel per_millisecond; K_rel=30; real V_rel micrometre_3; real Ca_rel(time) millimolar; when(time=time.min) Ca_rel=1.488; real u(time) dimensionless; when(time=time.min) u=2.35e-112; real v(time) dimensionless; when(time=time.min) v=1; real w(time) dimensionless; when(time=time.min) w=0.9992; real tau_u millisecond; real u_infinity(time) dimensionless; real tau_v(time) millisecond; real v_infinity(time) dimensionless; real tau_w(time) millisecond; real w_infinity(time) dimensionless; real i_tr(time) millimolar_per_millisecond; real tau_tr millisecond; tau_tr=180; real Ca_up(time) millimolar; when(time=time.min) Ca_up=1.488; real I_up_max millimolar_per_millisecond; I_up_max=0.005; real i_up(time) millimolar_per_millisecond; real K_up millimolar; K_up=0.00092; real i_up_leak(time) millimolar_per_millisecond; real Ca_up_max millimolar; Ca_up_max=15; real CMDN_max millimolar; CMDN_max=0.05; real TRPN_max millimolar; TRPN_max=0.07; real CSQN_max millimolar; CSQN_max=10; real Km_CMDN millimolar; Km_CMDN=0.00238; real Km_TRPN millimolar; Km_TRPN=0.0005; real Km_CSQN millimolar; Km_CSQN=0.8; real Ca_CMDN(time) millimolar; real Ca_TRPN(time) millimolar; real Ca_CSQN(time) millimolar; real V_cell micrometre_3; V_cell=20100; real V_i micrometre_3; real V_up micrometre_3; real B1(time) millimolar_per_millisecond; real B2(time) dimensionless; // // i_st=(if (((time>=stim_start) and (time<=stim_end)) and ((time-stim_start-floor((time-stim_start)/stim_period)*stim_period)<=stim_duration)) stim_amplitude else (0 picoA)); V:time=((-1)*(i_Na+i_K1+i_to+i_Kur+i_Kr+i_Ks+i_B_Na+i_B_Ca+i_NaK+i_CaP+i_NaCa+i_Ca_L+i_st)/Cm); // i_Na=(Cm*g_Na*m^3*h*j*(V-E_Na)); E_Na=(R*T/F*ln(Na_o/Na_i)); // alpha_m=(if (V=((-1)*(47.13 millivolt))) (3.2 per_millisecond) else (.32 per_millivolt_millisecond)*(V+(47.13 millivolt))/(1-exp((-1)*(.1 per_millivolt)*(V+(47.13 millivolt))))); beta_m=((.08 per_millisecond)*exp((-1)*V/(11 millivolt))); m_inf=(alpha_m/(alpha_m+beta_m)); tau_m=(1/(alpha_m+beta_m)); m:time=((m_inf-m)/tau_m); // alpha_h=(if (V<((-1)*(40 millivolt))) (.135 per_millisecond)*exp((V+(80 millivolt))/((-1)*(6.8 millivolt))) else (0 per_millisecond)); beta_h=(if (V<((-1)*(40 millivolt))) (3.56 per_millisecond)*exp((.079 per_millivolt)*V)+(3.1E5 per_millisecond)*exp((.35 per_millivolt)*V) else 1/((.13 millisecond)*(1+exp((V+(10.66 millivolt))/((-1)*(11.1 millivolt)))))); h_inf=(alpha_h/(alpha_h+beta_h)); tau_h=(1/(alpha_h+beta_h)); h:time=((h_inf-h)/tau_h); // alpha_j=(if (V<((-1)*(40 millivolt))) ((-1)*(127140.00000000001 per_millivolt_millisecond)*exp((.2444 per_millivolt)*V)-(3.4739999999999996E-5 per_millivolt_millisecond)*exp((-1)*(.04391 per_millivolt)*V))*(V+(37.78 millivolt))/(1+exp((.311 per_millivolt)*(V+(79.23 millivolt)))) else (0 per_millisecond)); beta_j=(if (V<((-1)*(40 millivolt))) (.1212 per_millisecond)*exp((-1)*(.01052 per_millivolt)*V)/(1+exp((-1)*(.1378 per_millivolt)*(V+(40.14 millivolt)))) else (.3 per_millisecond)*exp((-1)*(2.535E-7 per_millivolt)*V)/(1+exp((-1)*(.1 per_millivolt)*(V+(32 millivolt))))); j_inf=(alpha_j/(alpha_j+beta_j)); tau_j=(1/(alpha_j+beta_j)); j:time=((j_inf-j)/tau_j); // E_K=(R*T/F*ln(K_o/K_i)); i_K1=(Cm*g_K1*(V-E_K)/(1+exp((.07 per_millivolt)*(V+(80 millivolt))))); // i_to=(Cm*g_to*oa^3*oi*(V-E_K)); // alpha_oa=((.65 per_millisecond)*(exp((V-(-1)*(10 millivolt))/((-1)*(8.5 millivolt)))+exp((V-(-1)*(10 millivolt)-(40 millivolt))/((-1)*(59 millivolt))))^((-1)*1)); beta_oa=((.65 per_millisecond)*(2.5+exp((V-(-1)*(10 millivolt)+(72 millivolt))/(17 millivolt)))^((-1)*1)); tau_oa=((alpha_oa+beta_oa)^((-1)*1)/K_Q10); oa_infinity=((1+exp((V-(-1)*(10 millivolt)+(10.47 millivolt))/((-1)*(17.54 millivolt))))^((-1)*1)); oa:time=((oa_infinity-oa)/tau_oa); // alpha_oi=(((18.53 millisecond)+(1 millisecond)*exp((V-(-1)*(10 millivolt)+(103.7 millivolt))/(10.95 millivolt)))^((-1)*1)); beta_oi=(((35.56 millisecond)+(1 millisecond)*exp((V-(-1)*(10 millivolt)-(8.74 millivolt))/((-1)*(7.44 millivolt))))^((-1)*1)); tau_oi=((alpha_oi+beta_oi)^((-1)*1)/K_Q10); oi_infinity=((1+exp((V-(-1)*(10 millivolt)+(33.1 millivolt))/(5.3 millivolt)))^((-1)*1)); oi:time=((oi_infinity-oi)/tau_oi); // g_Kur=((.005 nanoS_per_picoF)+(.05 nanoS_per_picoF)/(1+exp((V-(15 millivolt))/((-1)*(13 millivolt))))); i_Kur=(Cm*g_Kur*ua^3*ui*(V-E_K)); // alpha_ua=((.65 per_millisecond)*(exp((V-(-1)*(10 millivolt))/((-1)*(8.5 millivolt)))+exp((V-(-1)*(10 millivolt)-(40 millivolt))/((-1)*(59 millivolt))))^((-1)*1)); beta_ua=((.65 per_millisecond)*(2.5+exp((V-(-1)*(10 millivolt)+(72 millivolt))/(17 millivolt)))^((-1)*1)); tau_ua=((alpha_ua+beta_ua)^((-1)*1)/K_Q10); ua_infinity=((1+exp((V-(-1)*(10 millivolt)+(20.3 millivolt))/((-1)*(9.6 millivolt))))^((-1)*1)); ua:time=((ua_infinity-ua)/tau_ua); // alpha_ui=(((21 millisecond)+(1 millisecond)*exp((V-(-1)*(10 millivolt)-(195 millivolt))/((-1)*(28 millivolt))))^((-1)*1)); beta_ui=((1 per_millisecond)/exp((V-(-1)*(10 millivolt)-(168 millivolt))/((-1)*(16 millivolt)))); tau_ui=((alpha_ui+beta_ui)^((-1)*1)/K_Q10); ui_infinity=((1+exp((V-(-1)*(10 millivolt)-(109.45 millivolt))/(27.48 millivolt)))^((-1)*1)); ui:time=((ui_infinity-ui)/tau_ui); // i_Kr=(Cm*g_Kr*xr*(V-E_K)/(1+exp((V+(15 millivolt))/(22.4 millivolt)))); // alpha_xr=(if (abs(V+(14.1 millivolt))<(1E-10 millivolt)) (.0015 per_millisecond) else (3E-4 per_millivolt_millisecond)*(V+(14.1 millivolt))/(1-exp((V+(14.1 millivolt))/((-1)*(5 millivolt))))); beta_xr=(if (abs(V-(3.3328 millivolt))<(1E-10 millivolt)) (3.7836118E-4 per_millisecond) else (7.3898E-5 per_millivolt_millisecond)*(V-(3.3328 millivolt))/(exp((V-(3.3328 millivolt))/(5.1237 millivolt))-1)); tau_xr=((alpha_xr+beta_xr)^((-1)*1)); xr_infinity=((1+exp((V+(14.1 millivolt))/((-1)*(6.5 millivolt))))^((-1)*1)); xr:time=((xr_infinity-xr)/tau_xr); // i_Ks=(Cm*g_Ks*xs^2*(V-E_K)); // alpha_xs=(if (abs(V-(19.9 millivolt))<(1E-10 millivolt)) (6.8E-4 per_millisecond) else (4E-5 per_millivolt_millisecond)*(V-(19.9 millivolt))/(1-exp((V-(19.9 millivolt))/((-1)*(17 millivolt))))); beta_xs=(if (abs(V-(19.9 millivolt))<(1E-10 millivolt)) (3.15E-4 per_millisecond) else (3.5E-5 per_millivolt_millisecond)*(V-(19.9 millivolt))/(exp((V-(19.9 millivolt))/(9 millivolt))-1)); tau_xs=(.5*(alpha_xs+beta_xs)^((-1)*1)); xs_infinity=((1+exp((V-(19.9 millivolt))/((-1)*(12.7 millivolt))))^((-1)*.5)); xs:time=((xs_infinity-xs)/tau_xs); // i_Ca_L=(Cm*g_Ca_L*d*f*f_Ca*(V-(65 millivolt))); // d_infinity=((1+exp((V+(10 millivolt))/((-1)*(8 millivolt))))^((-1)*1)); tau_d=(if (abs(V+(10 millivolt))<(1E-10 millivolt)) (4.579 millisecond)/(1+exp((V+(10 millivolt))/((-1)*(6.24 millivolt)))) else (1-exp((V+(10 millivolt))/((-1)*(6.24 millivolt))))/((.035 per_millivolt_millisecond)*(V+(10 millivolt))*(1+exp((V+(10 millivolt))/((-1)*(6.24 millivolt)))))); d:time=((d_infinity-d)/tau_d); // f_infinity=(exp((-1)*(V+(28 millivolt))/(6.9 millivolt))/(1+exp((-1)*(V+(28 millivolt))/(6.9 millivolt)))); tau_f=((9 millisecond)*(.0197*exp((-1)*(.0337 per_millivolt)^2*(V+(10 millivolt))^2)+.02)^((-1)*1)); f:time=((f_infinity-f)/tau_f); // f_Ca_infinity=((1+Ca_i/(3.5E-4 millimolar))^((-1)*1)); tau_f_Ca=(2 millisecond); f_Ca:time=((f_Ca_infinity-f_Ca)/tau_f_Ca); // sigma=(1/7*(exp(Na_o/(67.3 millimolar))-1)); f_NaK=((1+.1245*exp((-1)*.1*F*V/(R*T))+.0365*sigma*exp((-1)*F*V/(R*T)))^((-1)*1)); i_NaK=(Cm*i_NaK_max*f_NaK*1/(1+(Km_Na_i/Na_i)^1.5)*K_o/(K_o+Km_K_o)); // E_Ca=(R*T/(2*F)*ln(Ca_o/Ca_i)); i_B_Na=(Cm*g_B_Na*(V-E_Na)); i_B_Ca=(Cm*g_B_Ca*(V-E_Ca)); i_B_K=(Cm*g_B_K*(V-E_K)); // i_NaCa=(Cm*I_NaCa_max*(exp(gamma*F*V/(R*T))*Na_i^3*Ca_o-exp((gamma-1)*F*V/(R*T))*Na_o^3*Ca_i)/((K_mNa^3+Na_o^3)*(K_mCa+Ca_o)*(1+K_sat*exp((gamma-1)*V*F/(R*T))))); // i_CaP=(Cm*i_CaP_max*Ca_i/((5E-4 millimolar)+Ca_i)); // Fn=((1E3 millisecond_per_millimolar)*((1E-15 per_micrometre_3)*V_rel*i_rel-(1E-15 per_micrometre_3)/(2*F)*(.5*i_Ca_L-.2*i_NaCa))); i_rel=(K_rel*u^2*v*w*(Ca_rel-Ca_i)); // tau_u=(8 millisecond); u_infinity=((1+exp((-1)*(Fn-3.4175000000000003E-13)/1.367E-15))^((-1)*1)); u:time=((u_infinity-u)/tau_u); // tau_v=((1.91 millisecond)+(2.09 millisecond)*(1+exp((-1)*(Fn-3.4175000000000003E-13)/1.367E-15))^((-1)*1)); v_infinity=(1-(1+exp((-1)*(Fn-6.835E-14)/1.367E-15))^((-1)*1)); v:time=((v_infinity-v)/tau_v); // tau_w=(if (abs(V-(7.9 millivolt))<(1E-10 millivolt)) (6 millisecond)*.2/1.3 else (6 millisecond)*(1-exp((-1)*(V-(7.9 millivolt))/(5 millivolt)))/((1+.3*exp((-1)*(V-(7.9 millivolt))/(5 millivolt)))*(1 per_millivolt)*(V-(7.9 millivolt)))); w_infinity=(1-(1+exp((-1)*(V-(40 millivolt))/(17 millivolt)))^((-1)*1)); w:time=((w_infinity-w)/tau_w); // i_tr=((Ca_up-Ca_rel)/tau_tr); // i_up=(I_up_max/(1+K_up/Ca_i)); // i_up_leak=(I_up_max*Ca_up/Ca_up_max); // Ca_CMDN=(CMDN_max*Ca_i/(Ca_i+Km_CMDN)); Ca_TRPN=(TRPN_max*Ca_i/(Ca_i+Km_TRPN)); Ca_CSQN=(CSQN_max*Ca_rel/(Ca_rel+Km_CSQN)); // V_i=(V_cell*.68); V_rel=(.0048*V_cell); V_up=(.0552*V_cell); Na_i:time=(((-1)*3*i_NaK-(3*i_NaCa+i_B_Na+i_Na))/(V_i*F)); K_i:time=((2*i_NaK-(i_K1+i_to+i_Kur+i_Kr+i_Ks+i_B_K))/(V_i*F)); Ca_i:time=(B1/B2); B1=((2*i_NaCa-(i_CaP+i_Ca_L+i_B_Ca))/(2*V_i*F)+(V_up*(i_up_leak-i_up)+i_rel*V_rel)/V_i); B2=(1+TRPN_max*Km_TRPN/(Ca_i+Km_TRPN)^2+CMDN_max*Km_CMDN/(Ca_i+Km_CMDN)^2); Ca_up:time=(i_up-(i_up_leak+i_tr*V_rel/V_up)); Ca_rel:time=((i_tr-i_rel)*(1+CSQN_max*Km_CSQN/(Ca_rel+Km_CSQN)^2)^((-1)*1)); // }