У меня есть вопрос, касающийся реализации постоянно увеличивающегося вклада в программу, демонстрирующую мембранный потенциал нейрона в модели FitzHugh-Nagumo.
У меня есть структура модели, но я не знаюкак реализовать непрерывно увеличивающийся ввод.
import numpy as np
import matplotlib.pyplot as pl
*Fitzhugh-Nagumo model*
*responding to continuously increasing input*
# setup
#The initial value of the membrane potential (V0) should be -1.544
V0 = -1.544 # initial value for membrane potential (V)
#The initial value of the recovery variable (W0) should be -0.316
W0 = -0.316 # initial value for w
C = 1. # capacitance
b0 = 2. # intercept
b1 = 1.5 # slope
tau_w = 10. # time constant of w
#simlutate for 1 second
t_end = 1000. # end of simulation time (ms)
dt = 0.01 # length timestep (ms)
t_steps = int(t_end/dt)+1 # number of timesteps
T = np.linspace(0,t_end,t_steps) # array to keep track of time
v = np.linspace(-3.,3.,100) # a range of values for the membrane potential (V)
V = np.zeros(t_steps) # pre-allocate array to store membrane potential
W = np.zeros(t_steps) # pre-allocate array to store
V[0] = V0 # set first value to initial value (V)
W[0] = W0 # set first value to initial value (w)
#i need continuously increasing input
I = 2.5 # external input
# trajectory
for t in range(t_steps-1): # loop over time
dV = (V[t]-1./3.*np.power(V[t],3)-W[t]+I)/C # compute change of V
dW = (b0+b1*V[t]-W[t])/tau_w # compute change of w
V[t+1] = V[t]+dt*dV # update value of V
W[t+1] = W[t]+dt*dW # update value of w
# nullclines
w1 = I+ v-1./3.*np.power(v,3) # calculate w according to solution of dV/dt
w2 = b0+b1*v # calculate w according to solution of dw/dt
#plot membrane potential of the model over time
pl.figure(figsize=(20,20)) # create a figure of width and height equal to 10
pl.subplot(2,2,1) # fill first slot of 2-by-2 subplot (upper left)
pl.plot(v,w1,'r',label="dV/dt=0") # plot w1 as a function of v and label the line
pl.plot(v,w2,'b',label="dw/dt=0") # plot w2 as a function of v and label the line
pl.plot(V,W,'k',linewidth=2,label="trajectory") # plot trajectory of line defined by V and W over time, set line width, add label
pl.xlim(-3.,3.) # limit values of x-axis between -/+ 3
pl.ylim(-3.,5.) # limit values of y-axis between -3 and 5
pl.xlabel('V') # label the x-axis
pl.ylabel('w') # label the x-axis
pl.title('I=0') # give the subplot a title
pl.legend() # make the legend visible
pl.subplot(2,2,3) # fill third slot of 2-by-2 subplot (lower left)
pl.plot(T,V,'k') # plot membrane potential as a function of time
pl.xlim(0,t_end) # limit values of x-axis to the simulation time
pl.ylim(-3.,5.) # limit values of y-axis between -3 and 5
pl.xlabel('t') # label the x-axis
pl.ylabel('V') # label the y-axis
*input*
I = 2. # external input
*trajectory*
for t in range(t_steps-1): # loop over time
dV = (V[t]-1./3.*np.power(V[t],3)-W[t]+I)/C # compute change of u
dW = (b0+b1*V[t]-W[t])/tau_w # compute change of w
V[t+1] = V[t]+dt*dV # update value of u
W[t+1] = W[t]+dt*dW # update value of w
*nullclines*
w1 = I+ v-1./3.*np.power(v,3)
вычислить w в соответствии с решением du / dt w2 = b0 + b1 * v # вычислить w в соответствии с решением dw / dt
pl.subplot(2,2,2) # fill second slot of 2-by-2 *subplot (upper right)*
pl.plot(v,w1,'r',label="dV/dt=0") # plot w1 as a function of v and label the line
pl.plot(v,w2,'b',label="dw/dt=0") # plot w2 as a function of v and label the line
pl.plot(V,W,'k',linewidth=2,label="trajectory")
построить траекторию линии, определяемой V и W во времени, установить ширину линии, добавить метку pl.xlim (-3., 3.) # Предельные значения оси Xмежду - / + 3 pl.ylim (-3., 5.) # предельные значения оси y между -3 и 5 pl.xlabel ('V') # обозначить ось x pl.ylabel ('w')# обозначить ось x pl.title ('I = 2') # присвоить подзаголовку заголовок pl.legend () # сделать легенду видимой
pl.subplot(2,2,4) # fill third slot of 2-by-2 subplot (lower right)
pl.plot(T,V,'k') # plot membrane potential as a *function of time*
pl.xlim(0,t_end) # limit values of x-axis to the simulation time
pl.ylim(-3.,5.) # limit values of y-axis between -3 and 5
pl.xlabel('t') # label the x-axis
pl.ylabel('V') # label the y-axis
pl.show()
отобразить рисунок