#Nicoletti M, Loppini A, Chiodo L, Folli V, Ruocco G, et al. (2019) #Biophysical modeling of C. elegans neurons: Single ion currents and whole-cell dynamics of AWCon and RMD. #PLOS ONE 14(7): e0218738. https://doi.org/10.1371/journal.pone.0218738. PMID: 31260485 #=========================== RMD simulation ===================================================================== # potassium,calcium, leakge ad sodium reversal potentials par ek=-80,eca=60,eleak=-80,ena=30 #=============================================================================================================== # VOLTAGE-GATED POTASSIUM CHANNELS # ============================================================================================================== #==================== SHL-1 CHANNELS =========================================================================== #conductance in nS par gshal=2.48 par shalsfhit=18 # activation par vashal=11.2,kashal=14.1,kishal=8.3,vishal=-33.1 par ptmshal1=13.8,ptmshal2=-17.5165,ptmshal3=12.9213,ptmshal4=-3.7082,ptmshal5=6.4876,ptmshal6=1.8849 par cashal=0.1 minf_shal=1/(1+exp(-(v-vashal+shalsfhit)/(kashal))) tm_shal=(ptmshal1/(exp(-(v-ptmshal2)/ptmshal3)+exp((v-ptmshal4)/ptmshal5))+ptmshal6)*cashal #inactivation par pthfshal1=539.1584, pthfshal2=-28.1990, pthfshal3=4.9199, pthfshal4=27.2811 par pthsshal1=8422 pthsshal2=-37.7391 pthsshal3=6.3785 pthsshal4=118.8983,cshal=0.1 hinf_shal=1/(1+exp((v-vishal+shalsfhit)/(kishal))) ths_shal=(pthsshal1/(1+exp((v-pthsshal2)/pthsshal3))+pthsshal4)*cshal thf_shal=(pthfshal1/(1+exp((v-pthfshal2)/pthfshal3))+pthfshal4)*cshal m_shal'=(minf_shal-m_shal)/tm_shal hf_shal'=(hinf_shal-hf_shal)/thf_shal hs_shal'=(hinf_shal-hs_shal)/ths_shal init m_shal=0,hf_shal=1,hs_shal=1 Ishal=gshal*(m_shal^3)*(0.7*hf_shal+0.3*hs_shal)*(v-ek) #=================== SHK1 CHANNELS ============================================================================== #conductance par gshak=1.1 #activation par vashak=20.4,kashak=7.7 par ptmshak1=26.571450568169027 ,ptmshak2=-33.741611800716130,ptmshak3=15.757936311607475 par ptmshak4=15.364937728953288,ptmshak5=1.990037272604829 ,shiftV05=0 minf_shak=1/(1+exp(-(v-vashak+shiftV05)/kashak)) tm_shak=ptmshak1/(exp(-(v-(ptmshak2+shiftV05))/ptmshak4)+exp((v-(ptmshak2+shiftV05))/ptmshak3))+ptmshak5 #inactivation par kishak=5.8,vishak=-6.95 par pthshak=1400 hinf_shak=1/(1+exp((v-vishak+shiftV05)/kishak)) th_shak=pthshak m_shak'=(minf_shak-m_shak)/tm_shak h_shak'=(hinf_shak-h_shak)/th_shak init m_shak=0,h_shak=1 I_shak=gshak*m_shak*h_shak*(v-ek) #================ EGL-36 CHANNELS ============================================================================== #conductance par gegl36=1.3 #activation par va_egl36=63,ka_egl36=28.5,t1_egl36=355,t2_egl36=63,t3_egl36=13 par a1=0.31,a2=0.36,a3=0.39 min1_egl36=1/(1+exp(-(v-va_egl36)/(ka_egl36))) min2_egl36=1/(1+exp(-(v-va_egl36)/(ka_egl36))) min3_egl36=1/(1+exp(-(v-va_egl36)/(ka_egl36))) tm1_egl36=t1_egl36 tm2_egl36=t2_egl36 tm3_egl36=t3_egl36 m1_egl36'=(min1_egl36-m1_egl36)/tm1_egl36 m2_egl36'=(min2_egl36-m2_egl36)/tm2_egl36 m3_egl36'=(min3_egl36-m3_egl36)/tm3_egl36 init m1_egl36=0,m2_egl36=0,m3_egl36=0 I_egl36=gegl36*((a1*m1_egl36)+(a2*m2_egl36)+(a3*m3_egl36))*(v-ek) #======================= IRK CHANNELS ================================================================================= par gkir=0.2 par va_kir=-52, ka_kir=13 par p1tmkir=17.0752,p2tmkir=-17.8258, p3tmkir=20.3154, p4tmkir=-43.4414, p5tmkir=11.1691, p6tmkir=3.8329 minf_kir=1/(1+exp((v-va_kir+30)/ka_kir)) tm_kir=p1tmkir/(exp(-(v-p2tmkir)/p3tmkir)+exp((v-p4tmkir)/p5tmkir))+p6tmkir m_kir'=(minf_kir-m_kir)/tm_kir I_kir=gkir*m_kir*(v-ek) aux I_kir=I_kir #========================================================================================================================= # CALCIUM CHANNELS #========================================================================================================================= #=================== UNC-2 CHANNELS ======================================================================================= #conductance par gunc2=0.9 #activation par va_unc2=-12.17,ka_unc2=3.97,stm2=25 par p1tmunc2=1.4969,p2tmunc2=-8.1761,p3tmunc2=9.0753,p4tmunc2=15.3456,p5tmunc2=0.1029 par shiftmunc2=30,constmunc2=3 minf_unc2=1/(1+exp(-(v-va_unc2+stm2)/(ka_unc2))) tm_unc2=(p1tmunc2/(exp(-(v-p2tmunc2+shiftmunc2)/(p3tmunc2))+exp((v-p2tmunc2+shiftmunc2)/(p4tmunc2)))+p5tmunc2)*constmunc2 #inactivation par vi_unc2=-52.47,ki_unc2=5.6,sth2=25 par p1thunc2=83.8037, p2thunc2=52.8997,p3thunc2=3.4557,p4thunc2=72.0995,p5thunc2=23.9009,p6thunc2=3.5903 par consthunc2=1.7, shifthunc2=30 hinf_unc2=1/(1+exp((v-vi_unc2+sth2)/(ki_unc2))) th_unc2=(p1thunc2/(1+exp((v-p2thunc2+shifthunc2)/(p3thunc2)))+p4thunc2/(1+exp(-(v-p5thunc2+shifthunc2)/(p6thunc2))))*consthunc2 # equations m_unc2'=(minf_unc2-m_unc2)/tm_unc2 h_unc2'=(hinf_unc2-h_unc2)/th_unc2 init m_unc2=0,h_unc2=1 I_unc2=gunc2*m_unc2*h_unc2*(v-eca) #alpha and beta for slo1 ans slo2 channels alpha=minf_unc2/tm_unc2 beta=(1/tm_unc2)-alpha #================== EGL-19 CHANNELS =========================================================================================== #conductance par gegl19=0.99 #activation par va_egl19=5.6, ka_egl19=7.50,stm19=10,sth19=10 par pdg1=2.3359,pdg2=2.9324,pdg3=5.2357,pdg4=6.0,pdg5=1.8739,pdg6=1.3930,pdg7=30.0,stau19=10 minf_egl19=1/(1+exp(-(v-va_egl19+stm19)/ka_egl19)) tm_egl19=(pdg1+(pdg2*exp(-(v-pdg3+stau19)^2/(pdg4)^2))+(pdg5*exp(-(v-pdg6+stau19)^2/(pdg7)^2))) #inactivation par p1hegl19=1.4314,p2hegl19=24.8573,p3hegl19=11.9541,p4hegl19=0.1427,p5hegl19=5.9589,p6hegl19=-10.5428,p7hegl19=8.0552,p8hegl19=0.6038 par pds1=0.4,pds2=0.55,pds3=81.1179,pds4=-22.9723,pds5=5,pds6=43.0937,pds7=0.9,pds8=40.4885,pds9=28.7251,pds10=3.7125,pds11=0,shiftdps=10 hinf_egl19=((p1hegl19/(1+exp(-(v-p2hegl19+sth19)/p3hegl19))+p4hegl19)*(p5hegl19/(1+exp((v-p6hegl19+sth19)/p7hegl19))+p8hegl19)) ths_egl19=pds1*(((pds2*pds3)/(1+exp((v-pds4+shiftdps)/pds5)))+pds6+((pds7*pds8)/(1+exp((v-pds9+shiftdps)/pds10)))+pds11) # equations m_egl19'=(minf_egl19-m_egl19)/tm_egl19 hs_egl19'=(hinf_egl19-hs_egl19)/ths_egl19 init m_egl19=0,hs_egl19=1 I_egl19=gegl19*m_egl19*hs_egl19*(v-eca) # alpha and beta for slo1 and slo2 channels alpha1=minf_egl19/tm_egl19 beta1=(1/tm_egl19)-alpha1 #=============== CCA-1 channels ================================================================================================== #conductance par gcca1=3.1 par va_cca1=-42.65,ka_cca1=1.7,sscca1=15 par stmcca1=30,sthcca1=15,sshcca1=15,constmcca1=0.5 par fcca=1.4,f2cca1=1.15,f3ca=1.7,f4ca=1.1 par p1tmcca1=40,p2tmcca1=-62.5393,p3tmcca1=-12.4758,p4tmcca1=0.6947 minf_cca1=1/(1+exp(-(v-va_cca1+sscca1)/(ka_cca1*fcca))) tm_cca1=((p1tmcca1/(1+exp(-(v-p2tmcca1+stmcca1)/(p3tmcca1*f3ca))))+p4tmcca1)*constmcca1 par vi_cca1=-58,ki_cca1=7,consthcca1=0.08 par p1thcca1=280,p2thcca1=-60.7312,p3thcca1=8.5224,p4thcca1=19.7456 hinf_cca1=1/(1+exp((v-vi_cca1+sshcca1)/(ki_cca1*f2cca1))) th_cca1=((p1thcca1/(1+exp((v-p2thcca1+sthcca1)/(p3thcca1*f4ca))))+p4thcca1)*consthcca1 # equations m_cca1'=(minf_cca1-m_cca1)/tm_cca1 h_cca1'=(hinf_cca1-h_cca1)/th_cca1 init m_cca1=0,h_cca1=1 I_cca1=gcca1*m_cca1^2*h_cca1*(v-eca) #====================================================================================================================================== # CALCIUM-REGULATED POTASSIUM CHANNELS #====================================================================================================================================== #================== Ca nanodomain calculation (Only for SLO1 and SLO2)================================================================= # r: nanodomain radius # d: diffusion coefficient # F: Faraday constant C/mol # kb: # b: # gsc: conductance of calcium channels par r=13e-9,d=250e-12,F=96485,kb=500e6,b=30e-6,backgr=0.05,gsc=40e-12 cao_nano=(((abs(gsc*(v-eca)*1e-3)/(8*pi*r*d*F))*exp(-r/sqrt(d/(kb*b))))*1e6*1e-3 )+backgr cac_nano=backgr #================= SLO1 CHANNELS ==================================================================================================== # no clustering 1:1 stoichiometry with EGL-19 and UNC-2 channels # activation # (parameters from genetic algorithm) par wom=3.152961,wyx=0.012643,kyx=34.338784,nyx=0.000100 par wop=0.156217,wxy=-0.027527,kxy=55.726816,nxy=1.299198 kcm=wom*exp(-wyx*v)*(1/(1+((cac_nano)/kyx)^nyx)) kom=wom*exp(-wyx*v)*(1/(1+((cao_nano)/kyx)^nyx)) kop=wop*exp(-wxy*v)*(1/(1+(kxy/(cao_nano))^nxy)) kcp=0 # SLO1-UNC2 COMPLEX #conductance par gbk=0.3 minf_bk=(m_unc2*kop*(alpha+beta+kcm))/((kop+kom)*(kcm+alpha)+(beta*kcm)) tm_bkunc2=((alpha+beta+kcm)/((kop+kom)*(kcm+alpha)+(beta*kcm))) mbk'=(minf_bk-mbk)/tm_bkunc2 init mbk=0 I_bk=gbk*mbk*h_unc2*(v-ek) # SLO1-EGL19 COMPLEX #conductance par gslo1=0.3 #activation # same activation parameters of SLO1-UNC2 complex kcm2=wom*exp(-wyx*v)*(1/(1+((cac_nano)/kyx)^nyx)) kom2=wom*exp(-wyx*v)*(1/(1+((cao_nano)/kyx)^nyx)) kop2=wop*exp(-wxy*v)*(1/(1+(kxy/(cao_nano))^nxy)) kcp2=0 minf_slo1=(m_egl19*kop*(alpha1+beta1+kcm2))/((kop2+kom2)*(kcm2+alpha1)+(beta1*kcm2)) tm_slo1=((alpha1+beta1+kcm2)/((kop2+kom2)*(kcm2+alpha1)+(beta1*kcm2))) mslo1'=(minf_slo1-mslo1)/tm_slo1 init mslo1=0 I_slo1=gslo1*mslo1*hs_egl19*(v-ek) #================ SLO-2 CHANNELS ============================================================================================= # no clustering, 1:1 stoichiometry with EGL-19 and UNC-2 # activation # parameters from genetic algorithm par wom1=0.896395,wyx1=0.019405,kyx1=3294.553404,nyx1=0.000010 par wop1=0.026719,wxy1=-0.024123,kxy1=93.449423,nxy1=1.835067 kcm1=wom1*exp(-wyx1*v)*(1/(1+(cac_nano/kyx1)^nyx1)) kom1=wom1*exp(-wyx1*v)*(1/(1+((cao_nano)/kyx1)^nyx1)) kop1=wop1*exp(-wxy1*v)*(1/(1+(kxy1/(cao_nano))^nxy1)) kcp1=0 # SLO2- EGL19 COMPLEX #conductance par gbk2=0.3 minf_bk2=(m_egl19*kop1*(alpha1+beta1+kcm1))/((kop1+kom1)*(kcm1+alpha1)+(beta1*kcm1)) tm_bkegl19=((alpha1+beta1+kcm1)/((kop1+kom1)*(kcm1+alpha1)+(beta1*kcm1))) mbk2'=(minf_bk2-mbk2)/tm_bkegl19 init mbk2=0 I_bk2=gbk2*mbk2*hs_egl19*(v-ek) # SLO2-UNC2 COMPLEX #conductance par gslo2=0.3 kcm3=wom1*exp(-wyx1*v)*(1/(1+(cac_nano/kyx1)^nyx1)) kom3=wom1*exp(-wyx1*v)*(1/(1+((cao_nano)/kyx1)^nyx1)) kop3=wop1*exp(-wxy1*v)*(1/(1+(kxy1/(cao_nano))^nxy1)) kcp3=0 minf_slo2=(m_unc2*kop3*(alpha+beta+kcm3))/((kop3+kom3)*(kcm3+alpha)+(beta*kcm3)) tm_slo2=((alpha+beta+kcm3)/((kop3+kom3)*(kcm3+alpha)+(beta*kcm3))) mslo2'=(minf_slo2-mslo2)/tm_slo2 I_slo2=gslo2*mslo2*h_unc2*(v-ek) #================================= Intracellular calcium calculation ========================================================= # fd: Faraday constant in C/mmol par fd=96.485,fca=0.001,t_ca=50,rcell=1 I_ca=I_unc2+I_egl19+I_cca1 aux I_ca=I_ca # RMD total volume in um^3 (from Neuromorpho) par vol=5.65 alpha_ca=1/(2*vol*fd) J_ca1=fca*alpha_ca*I_ca aux J_ca1=J_ca1 backgr2=0.05e-3 rs=if(I_ca<0)then(-fca*(alpha_ca*I_ca)-((ca_intra1-backgr2)/t_ca))else((backgr2-ca_intra1)/t_ca) ca_intra1'=rs init ca_intra1=5e-05 #============================= KCNL CHANNELS ================================================================================== # Rinzel #conductance par gca=0.06 #activation par k_sk2=0.00033 minf_sk=ca_intra1/(k_sk2+ca_intra1) t_sk=6.3 #equations m_sk'=(minf_sk-m_sk)/t_sk init m_sk=0.13563 I_sk=gca*m_sk*(v-ek) #================================================================================================================================== # LEAKAGE CHANNELS #=================================================================================================================================== #================= LEAKAGE CURRENT ================================================================================================= par gleak=0.4 I_leak=gleak*(v-eleak) #================== Nca current ==================================================================================================== par gnca=0.05 I_nca=gnca*(v-ena) #=================================================================================================================================== # SIMULATIONS #==================================================================================================================================== # Total current I_tot=Ishal+I_bk+I_bk2+I_shak+I_cca1+I_egl19+I_unc2+I_sk+I_leak+I_egl36+I_nca+I_slo1+I_slo2+I_kir aux Itot=I_tot #===================== CURRENT CLAMP ================================================================================================ # IH: holding current # c: membrane capacitance pF par ih=0 par c=1.2 v'=(Istim-I_tot)/c init v=-70 #current clamp protocol Istim=if(t>ton)then(if(t<toff)then(iclamp)else(if(t>t2)then(if(t<t3)then(iclamp2)else(ih))else(ih)))else(ih) par iclamp=10,ton=310,toff=360,t2=410,iclamp2=-15,t3=430 aux prot=if(t>ton)then(if(t<toff)then(iclamp)else(if(t>t2)then(if(t<t3)then(iclamp2)else(ih))else(ih)))else(ih) ### XPP SETTINGS @ meth=stiff,trans=200,total=400,dt=0.01,maxstor=1000000,bound=100000000000000,atol=1e-8,tol=1e-8 @ xp=t,yp=v,xlo=200,xhi=1000,yhi=30,ylo=-100 done #===================== VOLTAGE CLAMP ================================================================================================= # voltage clamp protocol #par vk=-70 #v=if(t>ton)then(if(t<toff)then(vclamp)else(vk))else(vk) #par vclamp=1100,ton=1100,toff=2300 ### XPP SETTINGS #@ meth=stiff,trans=1000,total=2400,dt=0.01,maxstor=1000000,bound=100000000000000,atol=1e-8,tol=1e-8 #@ xp=t,yp=Itot,xlo=1000,xhi=2400,yhi=800,ylo=-100 #done