Neuronal Dynamics:Finite cable model 实验结果展示

一、演示效果和代码

通过设置输入的类型,在代码运行完之后会在本地产生一个显示输出内容的gif文件,和一个包含u-x和u-t图的png图片
当设置type=1的时候,相当于设定在t=0ms时在x=1输入一个冲激脉冲,会输出如下显示u-x随时间t变化的动态图:
请添加图片描述
当设置type=2的时候,相当于设定在t=0~1ms的时段内,在x=1持续输入的一个脉冲,会输出如下显示u-x随时间t变化的动态图:
请添加图片描述
当设置type=3的时候,相当于设定在t=0ms时在x=1输入一个冲激脉冲,同时在t=0~0.1ms的时段内,在x=4持续输入的一个脉冲,会输出如下显示u-x随时间t变化的动态图:
请添加图片描述

import numpy as np
from matplotlib import pyplot as plt
from matplotlib.pyplot import *
import matplotlib.animation
rc('font',size=16)
rc('font',family='SimHei')
plt.rcParams['axes.unicode_minus'] = False
class Cable():
    def __init__(self,L,all_time):
        self.L=L
        self.all_time=all_time
        self.dt=0.01
        self.dx=0.1
        self.I_list = [[0 for j in range(int(L/self.dx))] for i in range(int(all_time/self.dt))]
        self.x_list = [i*self.dx for i in range(int(L/self.dx))]
        self.t_list = [i*self.dt for i in range(int(all_time/self.dt))]
        self.u_list = [[0 for j in range(int(L/self.dx))] for i in range(int(all_time/self.dt))]
        self.type=0
    def change(self,type):
        self.type=type
        if type==0:
            self.I_list = [[0 for j in range(int(self.L/self.dx))] for i in range(int(self.all_time/self.dt))]
        elif type==1:
            self.I_list = [[(300 if (i*self.dt==0 and j==int(1/self.dx)) else 0) for j in range(int(self.L/self.dx))] for i in range(int(self.all_time/self.dt))]
        elif type==2:
            self.I_list = [[(10 if (i*self.dt<1 and j==int(1/self.dx)) else 0) for j in range(int(self.L/self.dx))] for i in range(int(self.all_time/self.dt))]
        elif type==3:
            self.I_list = [[(30 if (i!=0 and i*self.dt<0.1 and j==int(4/self.dx)) else (300 if (i*self.dt==0 and j==int(1/self.dx)) else 0)) for j in range(int(self.L/self.dx))] for i in range(int(self.all_time/self.dt))]
    def get_G(self,t,x):
        if t>0:
            return np.exp(-t-x*x/4/t)/np.sqrt(4*np.pi*t)
        return 0
    def get_G_limit(self,t0,x0,t,x):
        sum = 0
        for n in range(-10,10):
            sum+=self.get_G(t-t0,x-2*n*self.L-x0)+self.get_G(t-t0,x-2*n*self.L+x0)
        return sum
    def caculate_u(self):
        t0max = int(self.all_time/self.dt)
        x0max = int(self.L/self.dx)
        for x in range(x0max):
            for t in range(t0max):
                #通用计算过程
                for t0 in range(t):
                    for x0 in range(x0max):
                        if self.I_list[t0][x0]!=0:
                            self.u_list[t][x]+=self.get_G_limit(t0*self.dt,x0*self.dx,t*self.dt,x*self.dx)*self.I_list[t0][x0]
                #type1用下列代码计算更快
                # for t0 in range(0,1):
                #     self.u_list[t][x] += self.get_G_limit(t0*self.dt, 1, t*self.dt, x*self.dx) * 300
                self.u_list[t][x] *=self.dx*self.dt
def plot_line_pic(x,y,label='u-t变化图',ylim1=0,ylim2=0):
    plt.plot(x, y, '--r',linewidth=1.6, markeredgecolor='k', markersize=8)
    title(label,fontsize=20)
    if ylim1!=0 or ylim2!=0:
        ylim(ylim1,ylim2)
if __name__ == '__main__':
    model = Cable(L=5,all_time=1)
    # model.change(type=1)
    # model.change(type=2)
    model.change(type=3)
    model.caculate_u()
    subplot(121)
    if model.type==1:
        for i in range(3,40,6):
            plot_line_pic(model.x_list, model.u_list[i], label='u-x图;type='+str(model.type))
    else:
        for i in range(0,int(model.all_time/model.dt),int(model.all_time/model.dt/7)):
            plot_line_pic(model.x_list, model.u_list[i], label='u-x图;type='+str(model.type))
    subplot(122)
    for i in range(0,int(model.L/model.dx),(int((model.L/model.dx)/7))):
        plot_line_pic(model.t_list, [model.u_list[j][i] for j in range(len(model.t_list))], label='u-t图;type='+str(model.type))
    plt.savefig("Finite_cable_model_"+str(model.type)+".png")
    show()
    global fig
    fig, ax = plt.subplots()
    def update(i):
        plt.cla()
        plot_line_pic(model.x_list, model.u_list[int(i)], label='u-x图;type='+str(model.type)+';t='+str(round(i*model.dt,2)),ylim2=1)
    ani = matplotlib.animation.FuncAnimation(fig, update, frames=50, interval=100, repeat=True)
    ani.save("Finite_cable_model_"+str(model.type)+".gif", writer='imgemagick')

二、模型简述

该模型是对有限长度的链型神经结构的模拟:
由cable equation可以推出格林函数,即单位冲击函数对无穷长的链产生的影响为:
G ∞ ( t , x ) = Θ ( t ) 4 π t e x p ( − t − x 2 4 t ) G_\infty(t,x)=\frac{\Theta(t)}{\sqrt{4\pi t}}exp(-t-\frac{x^2}{4t}) G(t,x)=4πt Θ(t)exp(t4tx2)
考虑到限制函数:
∂ ∂ x u ( t , x = 0 ) = ∂ ∂ x u ( t , x = L ) = 0 \frac{\partial}{\partial x}u(t,x=0)=\frac{\partial}{\partial x}u(t,x=L)=0 xu(t,x=0)=xu(t,x=L)=0
由mirror charges技巧,单位冲击函数对有限长为L的链产生的影响为:
G 0 , L ( t 0 , x 0 ; t , x ) = ∑ n = − ∞ + ∞ G ∞ ( t − t 0 , x − 2 n L − x 0 ) + G ∞ ( t − t 0 , x − 2 n L + x 0 ) G_{0,L}(t_0,x_0;t,x)=\sum_{n=-\infty}^{+\infty}G_\infty(t-t_0,x-2nL-x_0)+G_\infty(t-t_0,x-2nL+x_0) G0,L(t0,x0;t,x)=n=+G(tt0,x2nLx0)+G(tt0,x2nL+x0)
结合输入单位时间和单位长度的外部电流输入iext(t,x)有:
u ( t , x ) = ∫ − ∞ t d t 0 ∫ 0 L d x 0 G 0 , L ( t 0 , x 0 ; t , x ) i e x t ( t 0 , x 0 ) u(t,x)=\int_{-\infty}^{t}dt_0\int_0^Ldx_0G_{0,L}(t_0,x_0;t,x)i_{ext}(t_0,x_0) u(t,x)=tdt00Ldx0G0,L(t0,x0;t,x)iext(t0,x0)

三、输出解释

3.1输入为t=0ms时位于x=1处的冲激函数

请添加图片描述
左图为0.03ms开始每隔0.06ms记录一次的u-x图结果
右图为在X=0与X=L之间均匀分布的七个位置的u随时间t的变化
在冲激函数的作用下,首先x=1处的幅度大幅上升,然后渐渐地下降,离x=1位置越远的地方摆动产生的幅度越小

3.2输入为t=0~1ms持续在x=1处输入脉冲

请添加图片描述
左图为在t=0与t=1之间均匀分布的七个时间点的u随时间x的变化
右图为在X=0与X=L之间均匀分布的七个位置的u随时间t的变化
在阶跃函数的作用下,各处的u值均持续上升,x=1处的幅度上升最大,离x=1位置越远的地方摆动产生的幅度越小

3.3输入为t=0ms位于x=1一个冲激脉冲,加上t=0~0.1ms的时段内,位于x=4持续输入的一个脉冲

请添加图片描述
可惜并没有出现和compartmental model一样的如下图的效果,推测是由于compartmental model考虑了分支的存在以及每一段的生物性质的差异,而我们这个模型中默认只有一条链
在这里插入图片描述

评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包

打赏作者

for-nothing

你的鼓励将是我创作的最大动力

¥1 ¥2 ¥4 ¥6 ¥10 ¥20
扫码支付:¥1
获取中
扫码支付

您的余额不足,请更换扫码支付或充值

打赏作者

实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值