CoCalc Public FilessimpleNeuronalCoupling.ipynbOpen with one click!
Authors: César Flores López, Roberto Garcia, Marco Herrera, Cesar Olvera, Fernando Perez, Jose Alberto Perez Benitez
Views : 164
Description: Tutorial to build simple networks of synaptically coupled neurons. The membrane potential dynamics are obtained using the thermodynamic transmembrane transport by Herrera-Valdez (PeerJ Preprints, 2017)
Compute Environment: Ubuntu 18.04 (Deprecated)

Neuronal dynamics in small networks: Simple synaptic coupling

Marco Arieli Herrera-Valdez1^1, Fernando Pérez Díaz2^2, Roberto García-Medina1^1, José Alberto Perez Benitez2^2

1^1 Facultad de Ciencias, Universidad Nacional Autónoma de México

2^2 ESIME-SEPI, Instituto Politécnico Nacional

Consider a general minimal model for transmembrane potential for a single cell where

(Share server only supports KaTeX; open in CoCalc to see this formula.)
Auxiliary functions
(Share server only supports KaTeX; open in CoCalc to see this formula.)

The amplitude of the synaptic current is described by a sum of an Ornstein-Uhlenbeck process representing non-specific inputs and a sum of synaptic currents that depend on the activity of the other cells. The parameters vxv_x, gxg_x, sxs_x, rxr_x represent, respectively, the half-activation potential, gating charge, de-activation bias, and basal rate for activation of the gate represented by the variable x{m,w}\quad x \in \left\{ m,w \right\}.

Simple network with plasticity

(Share server only supports KaTeX; open in CoCalc to see this formula.)
with auxiliary functions
(Share server only supports KaTeX; open in CoCalc to see this formula.)

\paragraph{Synapses between PINs} Let the connection matrix be w=(w11w1RwR1wRR)\begin{equation} w = \left( \begin{array}{ccc} w_{11}& \cdots & w_{1R} \\ \vdots& \ddots & \vdots \\ w_{R1}& \cdots & w_{RR} \end{array} \right) \end{equation} where

(Share server only supports KaTeX; open in CoCalc to see this formula.)
represents the connection weight from neuron jj to neuron ii. For simplification purposes, the activation of the presynaptic terminals in the iith neuron is assumed to depend on the transmembrane potential of the membrane potential of the iith neuron.

The excitatory postsynaptic current for the iith neuron is given by

(Share server only supports KaTeX; open in CoCalc to see this formula.)
where
(Share server only supports KaTeX; open in CoCalc to see this formula.)

In [ ]:
# Import necessary modules import scipy as sc import matplotlib.pylab as gr import time

Evolution of the system assuming that synapses activate as functions of time

In [ ]:
""" Change in membrane potential for one neuron. The function depends on the states of the system, represented by the vector U, time, and a dictionary containing parameters for the cell """ def NaTKaDNaKaAMPAGabaA(u,t,p): """ Notes: The voltages have been normalized y = v/vT and y* = v*/vT. Then dy/dt = dv/dt / vT The current amplitudes for the currents should be already normalized with respect to the membrane capacitance. A -> a = A/(vT Cm) """ y,w = u yyKaD = p['gain_Act_KaD']*(y-p['y_Half_Act_KaD']) wInf = expSigmoid_(yyKaD) wRate = p['rate_Act_KaD'] * expSum_(yyKaD,s=p['bias_Act_KaD']) mInf = expSigmoid_(p['gain_Act_NaT']*(y-p['y_Half_Act_NaT'])) jKaD= transmembraneFlux_(y,p["eta_KaD"]*p['y_Ka'], eta=p["eta_KaD"], rect=p["rect_KaD"], a=p["A_Single_KaD"]) jNaT= transmembraneFlux_(y,p["eta_NaT"]*p['y_Na'], eta=p["eta_NaT"], rect=p["rect_NaT"], a=p["A_Single_NaT"]) jNaKa= transmembraneFlux_(y,p['y_NaKa'], eta=p["eta_NaKa"], rect=p["rect_NaKa"], a=p["A_Single_NaKa"]) jAMPA= transmembraneFlux_(y,p['y_Ampa'], eta=p["eta_Ampa"], rect=p["rect_Ampa"], a=p["A_Single_Ampa"]) jGabaA= transmembraneFlux_(y,p['y_GabaA'], eta=p["eta_GabaA"], rect=p["rect_GabaA"], a=p["A_Single_GabaA"]) dy_Base = -p["N_NaT"]*mInf*(1-w)*jNaT - p["N_KaD"]*w*jKaD - p["N_NaKa"]* jNaKa dy_Syn = - p['N_Ampa']*p['p_Ampa'](t) * jAMPA - p['N_GabaA']*p['p_GabaA'](t) * jGabaA - p["A_LFP"](t) dw = (w**p['exp_Act_KaD']) * wRate * (wInf-w) return sc.array([dy_Base+dy_Syn,dw])

Funciones auxiliares

In [ ]:
# Function to calculate the transmembrane fluxes def transmembraneFlux_(y,y_r, eta, rect=0.5, a=1.0, degree="exp"): aa= (y_r- eta*y) A= a*eta if degree=="exp": i = A *((sc.exp((rect-1) * aa)-(sc.exp((rect* aa))))) elif degree=="3": i = A * aa * (1+ ((1-2* rect )/2)* aa + ((3* rect **2 -3* rect +1)/6)* (aa**2)) elif degree=="1": i = A * aa return i def expSum_(y, s=0.5): ee = sc.exp(y) return ee**(s-1)+ ee**(s) def expSigmoid_(y): ee = sc.exp(y) return ee/(1+ee) def wInverse_(w,vHalf,gCharge,vT): v= vHalf - (vT/gCharge) * sc.log( (1-w)/w) return v def expSub_(y, s=0.5): ee = sc.exp(y) return ee**(s-1) - ee**(s)

Solución numérica de ecuaciones

Change of variables to simplify the numerics

Let y=v/(CmvT)y= v/(C_m v_T) and y=v/(CmvT)y_* = v_*/(C_m v_T). Then

(Share server only supports KaTeX; open in CoCalc to see this formula.)
As a consequence, the amplitudes in the original system can also be normalized: NxaxAx=NxaxCmvT {N_x a_x} \mapsto A_x = \frac{N_x a_x}{C_m v_T}

In [ ]:
def normalizeVolts(parDict): """Normalize voltages with respect to the Boltzmann potential kT/q""" nDict={} vTCm= parDict['Cm'] * parDict['v_T'] nDict['vTCm']=vTCm for k,i in parDict.items(): if ((k.find('v_')==0)): #print(k) nDict["y_"+k[2:]]=parDict[k]/parDict['v_T'] parDict.update(nDict) return parDict # Normalization of potentials and amplitudes for the numerical solution of the system def normalizeAmps(parDict): """ Normalizes amplitudes for transmembrane currents with respect to the membrane capacitance. """ nDict={} vC= parDict['Cm'] * parDict['v_T'] nDict['vTCm']=vC for k,i in parDict.items(): if k.find('a_')==0: str0="A_" + k[2:] #if type(i)==float: nDict[str0]= i/vC #print(k,i); #print(str0, nDict[str0]) parDict.update(nDict) return parDict

Numerical method: Fourth order Runge-Kutta

In [ ]:
# ------------------------------------------ # Runge-Kutta numerics # ------------------------------------------ def RK4(f): return lambda y, t, dt: ( lambda dy1: ( lambda dy2: ( lambda dy3: ( lambda dy4: (dy1 + 2*dy2 + 2*dy3 + dy4)/6 )( dt * f( y + dy3, t + dt ) ) )( dt * f( y + dy2/2, t + dt/2 ) ) )( dt * f( y + dy1/2, t + dt/2 ) ) )( dt * f( y, t ) ) def solveRK4(p, f, ic, timeStep, nSteps): """ Example: u=solveRK4(p, f,ic,timeStep, nSteps) p is a dictionary containing the parameters for the system """ U=sc.zeros((nSteps, sc.prod(sc.shape(ic))),"float64") if len(p.keys()) * len(p.values()) >0: ff = lambda yy,tt: f(yy,tt,p) dU = RK4(ff) U[0]=ic for n in sc.arange(0,nSteps-1): changeU =dU( U[n], n*timeStep, timeStep) U[n+1] = U[n] + changeU return U.transpose()
In [ ]:
def calc2DNrnOrbit(p,graphs=1,fs=16, maxT=20, show_wInv=0, plotCurrents=0, plotTimeCurrents=0, figName=''): """ Examples: q=calc2DNrnOrbit(p,graphs=0) q=calc2DNrnOrbit(p,graphs=1) """ orbit=solveRK4(p=p, f=p["rhs"],ic=p["ic"], timeStep=p["timeStep"], nSteps=p["nSteps"]) y = orbit[0]; v = y*p["v_T"] w = orbit[1] dvdt=sc.zeros(len(v)); dvdt[1:]= (v[1:]-v[:-1])/p["timeStep"] dwdt=sc.zeros(len(v)); dwdt[1:]= (w[1:]-w[:-1])/p["timeStep"] #wInv= wInverse(w,vHalf=p["v_Half_Act_KaD"],vT=p["v_T"],gCharge=p["gain_Act_KaD"]) print("max dv/dt=%g"%(dvdt.max())) #print("orbit: ", orbit) yyKaD = p['gain_Act_KaD']*(y-p['y_Half_Act_KaD']) wInf = expSigmoid_(yyKaD) wRate = p['rate_Act_KaD'] * expSum_(yyKaD,s=p['bias_Act_KaD']) mInf = expSigmoid_(p['gain_Act_NaT']*(y-p['y_Half_Act_NaT'])) jKaD= transmembraneFlux_(y,p['y_Ka'], eta=p["eta_KaD"], rect=p["rect_KaD"], a=p["a_Single_KaD"],degree="exp") jNaT= transmembraneFlux_(y,-p['y_Na'], eta=p["eta_NaT"], rect=p["rect_NaT"], a=p["a_Single_NaT"],degree="exp") jNaKa= transmembraneFlux_(y,p['y_NaKa'], eta=p["eta_NaKa"], rect=p["rect_NaKa"], a=p["a_Single_NaKa"],degree="exp") iNaT = p["N_NaT"]*(1-w)*mInf *jNaT iKaD = p["N_KaD"]*w* jKaD iNaKa = p["N_NaKa"]*jNaKa iTot= (iNaT+ iKaD + iNaKa) upInds = sc.where(dvdt>0)[0] dnInds = sc.where(dvdt<0)[0] rcNaT = iNaT.sum()/iTot.sum() rcKd = iKaD.sum()/iTot.sum() effNaT = iNaT[upInds].sum()/iTot[upInds].sum() effKd = iKaD[dnInds].sum()/iTot[dnInds].sum() print("Na efficiency = %g"%((effNaT))) print("K efficiency = %g"%((effKd))) #jAMPA = p["Cm"]*p['AMPA'](timeSamp) * nExpSub(y- p['y_r_AMPA'], s=p['bias_AMPA']) #jGabaA = p["Cm"]*p['GabaA'](timeSamp) * nExpSub(y- p['y_r_GabaA'], s=p['bias_GabaA']) p['timeSamples'] =sc.arange(0,p['timeMax'],p['timeStep']) q={"timeSamples":p["timeSamples"], "v":v, "w":w, "dvdt":dvdt, "dwdt":dwdt, #"wInv":wInv, "mInf":mInf, "wInf":wInf, "iNaT":iNaT, "iKaD":iKaD, "iNaKa":iNaKa, "iTot":iTot, #"vNull":vNull, #"jAMPA":jAMPA, "jGabaA":jGabaA } if graphs>0: maxC = sc.array([iNaT.max(),iKaD.max()]).max() f1=gr.figure(figsize=(17,5)); q["fig1"]=f1; if plotCurrents>0: f2=gr.figure(figsize=(17,5));q["fig3"]=f3; q["fig2"]=f2 r=2; c=2; ax=list() for n in range(r*c): ax.append(f2.add_subplot(r,c,n+1)) gr.ioff() # Generar los planos para las graficas axV=f1.add_subplot(1,1,1) # Asignar a cada plano una grafica #xx=wInv*p["rate_Act_KaD"] axV.plot(p['timeSamples'], v, 'k', lw=1, alpha=1.0, label=r'$(t,v)$'); ax[0].plot(dvdt, v, 'k', lw=1, alpha=1.0, label=r'$(\partial_t v,v)$'); ax[0].plot([dvdt.min(),dvdt.max()], [0,0], 'k:'); ax[1].plot(w, v, 'k', lw=1, alpha=1.0, label=r'$(w,v)$'); if show_wInv>0: #axV.plot(p['timeSamples'], wInv, 'k', lw=3, alpha=0.35, label=r'$(t,w^{-1}(v))$'); ax[1].plot(w, wInv, 'k', lw=3, alpha=0.35, label=r'$(w,w^{-1}(v))$'); ax[0].plot(dvdt, wInv, 'k', lw=3, alpha=0.3, label=r'$(\partial_t v,w^{-1}(v))$'); if plotTimeCurrents>0: ax.plot(p['timeSamples'], iKaD, 'r', lw=2, alpha=1.0, label=r'$(t,i_{Kd})$'); ax.plot(p['timeSamples'], iNaT, 'g', lw=2, alpha=1.0, label=r'$(t,i_{NaT})$'); ax.plot(p['timeSamples'], iNaKa, 'b', lw=2, alpha=1.0, label=r'$(t,i_{NaK})$'); ax.plot(p['timeSamples'],iTot, 'k', lw=1, alpha=1.0, label=r'$(t,i_{M})$'); ax.plot([0,p['timeMax']], [0,0], 'k:'); ax.legend(ncol=1,loc="upper right",fontsize=fs) ax.set_xlim(p["timeSamples"].min(),maxT) ax.set_ylabel('nA',fontsize=fs); axC.set_xlabel('ms',fontsize=fs) f3.text(xStart,0.975,'B',fontsize=fs) ax[2].plot(dvdt, iKaD, 'r', lw=2, alpha=1.0, label=r'$(\partial_t v,i_{Kd})$'); ax[2].plot(dvdt, iNaT, 'g', lw=2, alpha=1.0, label=r'$(\partial_t v,i_{NaT})$'); ax[2].plot(dvdt, iNaKa, 'b', lw=2, alpha=1.0, label=r'$(\partial_t v,i_{NaK})$'); ax[2].plot(dvdt, iTot, 'k', lw=1, alpha=1.0, label=r'$(\partial_t v,i_{M})$'); ax[2].plot([dvdt.min(),dvdt.max()], [0,0], 'k:'); ax[3].plot(w, iKaD, 'r', lw=2, alpha=1.0, label=r'$(\partial_t w,i_{Kd})$'); ax[3].plot(w, iNaT, 'g', lw=2, alpha=1.0, label=r'$(\partial_t w,i_{NaT})$'); ax[3].plot(w, iNaKa, 'b', lw=2, alpha=1.0, label=r'$(\partial_t w,i_{NaK})$'); ax[3].plot(w, iTot, 'k', lw=1, alpha=1.0, label=r'$(\partial_t w,i_{M})$'); ax[3].plot([w.min(),w.max()], [0,0], 'k:'); axV.legend(ncol=1,loc="upper right",fontsize=fs) ax[0].legend(ncol=1,loc="lower right",fontsize=fs) ax[1].legend(ncol=1,loc="lower right",fontsize=fs) ax[2].legend(ncol=1,loc="upper right",fontsize=fs) ax[3].legend(ncol=1,loc="upper left",fontsize=fs) ax[2].set_ylim(-maxC,maxC) ax[3].set_ylim(-maxC,maxC) axV.set_ylabel('$v$ (mV)',fontsize=fs); axV.set_xlabel('ms',fontsize=fs); axV.set_xlim(p["timeSamples"].min(),maxT); ax[0].set_ylabel('$v$ (mV)',fontsize=fs); ax[0].set_xlabel('V/s',fontsize=fs); ax[1].set_ylabel('$v$ (mV)',fontsize=fs); ax[1].set_xlabel('$w$',fontsize=fs); ax[2].set_ylabel('pA',fontsize=fs); ax[2].set_xlabel('V/s',fontsize=fs) ax[3].set_ylabel('pA',fontsize=fs); ax[3].set_xlabel('$\partial_t w$ (1/s)',fontsize=fs) f1.subplots_adjust(left=0.05, bottom=0.05, right=0.95, top=0.95, wspace=0.2, hspace=0.13) f2.subplots_adjust(left=0.05, bottom=0.05, right=0.95, top=0.95, wspace=0.2, hspace=0.13) xStart=0.01; xDiv=0.485; yDiv=0.485 f1.text(xStart,0.975,'A',fontsize=fs) f2.text(xStart,0.975,'C',fontsize=fs) f2.text(xStart+xDiv,0.975,'D',fontsize=fs) f2.text(xStart,yDiv,'E',fontsize=fs) f2.text(xStart+xDiv,yDiv,'F',fontsize=fs) gr.ion(); gr.draw(); if len(figName)>0: n1=figName[:-4]+'_p1.png' n2=figName[:-4]+'_p2.png' f1.savefig(n1) f2.savefig(n2) print("Saved figure to files %s and %s"%(n1,n2)) if plotCurrents: n3=figName[:-4]+'_currents.png' f3.savefig(n3) print("and the currents to %s and %s"%(n3)) return q
In [ ]:

Parameters

In [ ]:
# Dictionaries with parameters and functions pars={'Cm':20.0, 'v_T':26.7, 'gain_Act_NaT':5.0, 'v_Half_Act_NaT': -17.0, 'gain_Act_KaD':3.5, 'v_Half_Act_KaD': 0.0, 'rate_Act_KaD': 2.0, 'bias_Act_KaD': 0.3,'exp_Act_KaD':1.0, 'v_Ka':-90.0, 'v_Na':60.0, 'v_ATP':-430.0, 'v_Ampa':0.0, 'v_GabaA':-70.0, 'a_Single_NaT':1.2, 'a_Single_KaD':1.0,'a_Single_NaKa':1.0, 'a_Single_Ampa':1.0,'a_Single_GabaA':1.0, 'N_NaT':1400.0, 'N_KaD':4500.0, 'N_NaKa':200.0, 'N_Ampa':0.0, 'N_GabaA':0.0, 'eta_NaT':-1.0, 'eta_KaD':1.0, 'eta_NaKa':1.0, 'eta_Ampa':-1.0, 'eta_GabaA':1.0, 'rect_KaD':0.5,'rect_NaT':0.5,'rect_NaKa':0.5,'rect_Ampa':0.5,'rect_GabaA':0.5, } pars["v_NaKa"]= 3*pars["v_Na"] - 2*pars["v_Ka"] + pars["v_ATP"] funcs ={"LFP": lambda t: 0.0} descriptions ={'Assumptions':"Potassium KD channel activation describes NaT channel inactivation.", 'gain_Act_NaT': "Ganancia en el estado estable de activación de canales de sodio", 'gain_Act_KaD': "Ganancia en el estado estable de activación de canales de potasio", 'v_Half_Act_NaT': "Voltaje de activación media de canales de sodio", 'v_Half_Act_KaD': "Voltaje de activación media de canales de potasio", 'v_ATP': "Voltaje de reversa de la ATPasa de Na-K"}
Calculation of the solution

Initial conditions and number of steps

In [ ]:
pars["nCells"]=1 v0=-30.0/pars["v_T"] y0= 0.001 pars["ic"]= sc.array([v0,y0]) pars["timeMax"] = 40.0 pars["timeStep"]=1/70.0; pars["nSteps"]=10000; pars["nSteps"]= sc.int32(pars["timeMax"]/pars["timeStep"]) pars["timeSamples"]=sc.arange(0,pars["timeMax"],pars["timeStep"]) pars['N_Ampa']=10.0; pars["p_Ampa"]= lambda t: sc.int32((t>20)&(t<22)) + sc.int32((t>30)&(t<32)) pars["p_GabaA"]= lambda t: 0.0 pars["A_LFP"] = lambda t: 0.0 pars["rhs"]=NaTKaDNaKaAMPAGabaA pars= normalizeAmps(pars) pars= normalizeVolts(pars)

The solution to the system is composed of two vectors that correspond to the two variables of the system, yy and ww. Then y can be transformed back into vv by setting v=yvTv = y\cdot v_T.

In [ ]:
orbit=solveRK4(p=pars, f=pars["rhs"],ic=pars["ic"], timeStep=pars["timeStep"], nSteps=pars["nSteps"]) y = orbit[0]; v = y*pars["v_T"] w = orbit[1] dvdt=sc.zeros(len(v)); dvdt[1:]= (v[1:]-v[:-1])/pars["timeStep"] f1=gr.figure(figsize=(15,5)) ax_tv=f1.add_subplot(121); ax_dvdt=f1.add_subplot(122) ax_tv.plot(pars["timeSamples"],v,label=r'$(t,v)$') ax_dvdt.plot(dvdt,v,label=r'$(\partial_t v,v)$') ax_tv.legend(); ax_dvdt.legend() #q=calc2DNrnOrbit(pars,graphs=1,fs=16, maxT=40, show_wInv=0, plotCurrents=0, figName='')

Network dynamics with simple synaptic coupling

The activation of synapses depends only on the voltage of the presynaptic neuron. Assume for simplicity that there are N=2N=2 neurons in the network. Construct two matrices for the connections such that the first neuron excites the second, and the second neuron excites the first: CAmpa=[0100],CGabaA=[0010]\begin{equation} C_{Ampa} =\left[ \begin{array}{cc} 0&1\\0&0 \end{array} \right], \quad C_{GabaA} =\left[ \begin{array}{cc} 0&0\\1&0 \end{array} \right] \end{equation}

In [ ]:
pars["gain_Act_AMPA"] = 5.0; pars["v_Half_Act_AMPA"] = -20.0 pars["gain_Act_GabaA"] = 5.0; pars["v_Half_Act_GabaA"] = -20.0 pars=normalizeVolts(pars)
In [ ]:
# Connection matrices pars["conn_Ampa"]= sc.matrix([[0,1],[0,0]]) print(pars["conn_Ampa"]) print(sc.array([1,2])) sc.dot(pars["conn_Ampa"],sc.array([1,2]))
In [ ]:
""" Change in membrane potential for one neuron. The function depends on the states of the system, represented by the vector U, time, and a dictionary containing parameters for the cell """ def NaTKaDNaKaAMPAGabaA_Net(u,t,p): """ Notes: The voltages have been normalized y = v/vT and y* = v*/vT. Then dy/dt = dv/dt / vT The current amplitudes for the currents should be already normalized with respect to the membrane capacitance. A -> a = A/(vT Cm) """ y= u[0:pars["nCells"]] w= u[pars["nCells"]: 2*pars["nCells"]] yyKaD = p['gain_Act_KaD']*(y-p['y_Half_Act_KaD']) wInf = expSigmoid_(yyKaD) wRate = p['rate_Act_KaD'] * expSum_(yyKaD,s=p['bias_Act_KaD']) mInf = expSigmoid_(p['gain_Act_NaT']*(y-p['y_Half_Act_NaT'])) jKaD= transmembraneFlux_(y,p["eta_KaD"]*p['y_Ka'], eta=p["eta_KaD"], rect=p["rect_KaD"], a=p["A_Single_KaD"]) jNaT= transmembraneFlux_(y,p["eta_NaT"]*p['y_Na'], eta=p["eta_NaT"], rect=p["rect_NaT"], a=p["A_Single_NaT"]) jNaKa= transmembraneFlux_(y,p['y_NaKa'], eta=p["eta_NaKa"], rect=p["rect_NaKa"], a=p["A_Single_NaKa"]) jAMPA= transmembraneFlux_(y,p['y_Ampa'], eta=p["eta_Ampa"], rect=p["rect_Ampa"], a=p["A_Single_Ampa"]) jGabaA= transmembraneFlux_(y,p['y_GabaA'], eta=p["eta_GabaA"], rect=p["rect_GabaA"], a=p["A_Single_GabaA"]) dw = (w**p['exp_Act_KaD']) * wRate * (wInf-w) pAmpa_ss = sc.dot(pars["conn_Ampa"], expSigmoid_(p['gain_Act_Ampa']*(y-p['y_Half_Act_Ampa']))) pGabaA_ss = sc.dot(pars["conn_GabaA"], expSigmoid_(p['gain_Act_GabaA']*(y-p['y_Half_Act_GabaA']))) pGabaA_ss = expSigmoid_(p['gain_Act_GabaA']*(y-p['y_Half_Act_GabaA'])) dy_Base = -p["N_NaT"]*mInf*(1-w)*jNaT - p["N_KaD"]*w*jKaD - p["N_NaKa"]* jNaKa dy_Syn = - p['N_Ampa']*p['p_Ampa'](t) * jAMPA - p['N_GabaA']*p['p_GabaA'](t) * jGabaA - p["A_LFP"](t) return sc.array([dy_Base+dy_Syn,dw])
In [ ]:
pars["nCells"]=2 v0=-sc.ones(pars["nCells"])*30.0/pars["v_T"] y0= sc.ones(pars["nCells"]) pars["ic"]= sc.array([v0,y0]) print(pars["ic"])
In [2]:
pars["nCells"]=2 v0=-sc.ones(pars["nCells"])*30.0/pars["v_T"] y0= sc.ones(pars["nCells"]) pars["ic"]= sc.array([v0,y0]) pars["timeMax"] = 40.0 pars["timeStep"]=1/70.0; pars["nSteps"]=10000; pars["nSteps"]= sc.int32(pars["timeMax"]/pars["timeStep"]) pars["rhs"]=NaTKaDNaKaAMPAGabaA pars= normalizeAmps(pars) pars= normalizeVolts(pars)
--------------------------------------------------------------------------- NameError Traceback (most recent call last) <ipython-input-2-731127457a27> in <module>() ----> 1 pars["nCells"]=2 2 v0=-sc.ones(pars["nCells"])*30.0/pars["v_T"] 3 y0= sc.ones(pars["nCells"]) 4 pars["ic"]= sc.array([v0,y0]) 5 pars["timeMax"] = 40.0 NameError: name 'pars' is not defined
In [ ]:
In [ ]: