Effect of polysynaptic facilitaiton between piriform-hippocampal network stages (Trieu et al 2015)

 Download zip file 
Help downloading and running models
Accession:181032
This is a model of a multistage network with stages representing regions and synaptic contacts from the olfactory cortex to region CA1 of the hippocampus in Brian2 spiking neural network simulator (Trieu et al 2015). It is primarily designed to assess how synaptic facilitation at multiple stages in response to theta firing changes the output of the network. Further developments will be posted at: github.com/cdcox/multistage_network This model was prepared by Conor D Cox, University of California, Irvine For questions please contact Conor at cdcox1@gmail.com
Reference:
1 . Trieu BH, Kramár EA, Cox CD, Jia Y, Wang W, Gall CM, Lynch G (2015) Pronounced differences in signal processing and synaptic plasticity between piriform-hippocampal network stages: a prominent role for adenosine. J Physiol 593:2889-907 [PubMed]
Citations  Citation Browser
Model Information (Click on a link to find other models with that property)
Model Type: Realistic Network;
Brain Region(s)/Organism: Hippocampus; Dentate gyrus; Olfactory cortex; Generic; Entorhinal cortex;
Cell Type(s):
Channel(s):
Gap Junctions:
Receptor(s):
Gene(s):
Transmitter(s):
Simulation Environment: Brian; Python;
Model Concept(s): Synaptic Integration; Olfaction;
Implementer(s): Cox, Conor [cdcox1 at gmail.com];
/
TrieuEtAl2015
README.txt
5hzresponsedata.csv
multistage_network_model_time.py
                            
# -*- coding: utf-8 -*-
"""
Created on Mon Dec 08 15:24:48 2014

@author: Conor D. Cox For paper Pronounced differences in signal processing 
and synaptic plasticity between piriform-hippocampal network stages:
A prominent role for adenosine Figure 9 last panel
We are improving the code and updates will be posted to github.com/cdcox/multistage_network
This is a multi-stage network model to demonstrate the changes a signal goes through
due to regional facilitation in response to theta firing.
For questions please contact Conor D Cox at cdcox1@gmail.com
"""

from brian2 import *
import numpy as np

def smooth(y, box_pts):
    '''used to smooth regions of each graph against random variations'''
    box = np.ones(box_pts)/box_pts
    y_smooth = convolve(y, box, mode='same')
    return y_smooth
    
def string_assemble(data):
    '''takes input graphs, smooths the, then finds a derivative point far enough away
    for derivative to be 'meaningful' then builds piecewise function based off critical points'''
    data=[float(x) for x in data if x!='']
    data=data[:100]
    data=smooth(data,10)
    data=data[:-10]
    derivative=data[7:]-data[:-7]
    edges=np.where((derivative>0)[1:]-(derivative>0)[:-1])
    output_str='0'
    first_edge=0
    for NN,number in enumerate(edges[0]):
        fit=polyfit(np.arange(first_edge,number)*200,data[first_edge:number],1)
        output_str=output_str+r'+(('+str(fit[0])+r')*t/ms+'+str(fit[1])+r')'+r'*(t/ms>'+str(first_edge*200)+r')*(t/ms<='+str(number*200)+r')'
        first_edge=number
    fit=polyfit(np.arange(first_edge,len(data))*200,data[first_edge:len(data)],1)
    output_str=output_str+r'+(('+str(fit[0])+r')*t/ms+'+str(fit[1])+r')'+r'*(t/ms>'+str(first_edge*200)+r')'
    return(output_str)
    
def build_model_string(equation,firing_strength):
    '''Takes piecewise functions and sets them in a format appropriate for Brian2'''
    string_out=r'scaler='+equation+'''
    x=x+1 
    ge+='''+str(firing_strength)+'''*scaler/100.*mV'''
    return(string_out)

'''Opens book of initial graphs and clean them up for use in Brain2'''    
book = np.genfromtxt(r'5hzresponsedata.csv',delimiter=',',skiprows=3)

#############################################################################
LOT_data=np.array(book[0:96,1])
LOT_str=string_assemble(LOT_data)
 #############################################################################   
AP_data=np.array(book[0:96,2])
AP_str=string_assemble(AP_data)
#############################################################################
PE_data=np.array(book[0:96,3])
PE_str=string_assemble(PE_data)
#############################################################################
LPP_data=np.array(book[0:96,4])
LPP_str=string_assemble(LPP_data)
#############################################################################
DG_data=np.array(book[0:96,5])
DG_str=string_assemble(DG_data)
#############################################################################
CA3_recurrent_data=np.array(book[0:96,6])
CA3_str=string_assemble(CA3_recurrent_data)
#############################################################################
CA1_data=np.array(book[0:96,6])
CA1_str=string_assemble(CA1_data)
#############################################################################
'''Define connectivity matrix based on IO curve values'''
names=['LOT','ASSN','PE','LPP','MF','CA1','CA3',]
x=[1.8,1.8,1.8,1.8,-1,.2,.2]
x=[1-(stuff)/3.1/1.7 for stuff in x]
fix=.15
con_list=np.array(x)*fix
'''baseline equation of excitatory and inhibitory decays'''
eqs = '''
dv/dt = ((ge+gi)/(20.*ms)-(v+60.*mV)/(20.*ms)) : volt
dge/dt = -ge/(5.*ms) : volt
dgi/dt = -gi/(10.*ms) : volt
'''

'''Book keeping to set up model'''
generalized_variable_model='''x:1
                              scaler:1'''
strength=3.2 #excitatory strength
'''Setting up olfactory bulb spiking'''
start=0
indices=[]
times=[]

for i in range(200):
    for j in range(200):
        indices.append(j)
        times.append(start)
    start+=200
indices = array(indices)
times = array(times - 0.5*defaultclock.dt/ms)*ms #moving spiking times to fix Brian2 bug
Olfactory_bulb = SpikeGeneratorGroup(200, indices, times)

'''Each region is set up with its number of cells, equation of firing strength change over time
and connectivity and after hyper-polarization'''
Piriform=NeuronGroup(200, eqs, threshold='v>-50*mV', reset='v=-60*mV')
Piriform.v = -60*mV
model_string=build_model_string(LOT_str,strength)
LOT = Synapses(Olfactory_bulb, Piriform, model=generalized_variable_model,pre=model_string)
LOT.connect(True, p=con_list[0])
##############################################################################
model_string=build_model_string(AP_str,strength)
CP_P=Synapses(Piriform, Piriform, model=generalized_variable_model,pre=model_string)                                            
CP_P.connect(True,  p=con_list[1])
CP_P.delay=1*ms  
##########################################################################################
CP_ah=Synapses(Piriform,Piriform,pre=r'gi-=4*mV')
for i in range(200):
    CP_ah.connect(i,i)
#############################################################################################
Entorhinal=NeuronGroup(200, eqs, threshold='v>-50*mV', reset='v=-60*mV') 
Entorhinal.v = -60*mV
model_string=build_model_string(PE_str,strength)
CP_e= Synapses(Piriform, Entorhinal,model=generalized_variable_model,pre=model_string)
CP_e.connect(True, p=con_list[2])
#############################################################################
CE_ah=Synapses(Entorhinal,Entorhinal,pre=r'gi-=4*mV')
for i in range(200):
    CE_ah.connect(i,i)
############################################################################
Dentate=NeuronGroup(200, eqs, threshold='v>-50*mV', reset='v=-60*mV')
Dentate.v = -60*mV
model_string=build_model_string(LPP_str,strength)
LPP = Synapses(Entorhinal, Dentate,model=generalized_variable_model,pre=model_string)
LPP.connect(True, p=con_list[3])
##############################################################################
CD_ah=Synapses(Dentate,Dentate,pre=r'gi-=4*mV')
for i in range(200):
    CD_ah.connect(i,i)
############################################################################
CA3=NeuronGroup(400, eqs, threshold='v>-50*mV', reset='v=-60*mV')
CA3.v = -60*mV
model_string=build_model_string(DG_str,10)
MF = Synapses(Dentate, CA3,model=generalized_variable_model,pre=model_string)
MF.connect(True, p=con_list[4])
###########################################################################
'''simulation of simplified CA3 inhibition'''
C3_ah=Synapses(CA3,CA3,pre=r'gi-=4*mV')
for i in range(400):
    C3_ah.connect(i,i)
    
C3_recurrents=Synapses(CA3,CA3,pre=r'gi-=2*mV')
C3_recurrents.connect(True,p=.15)
C3_recurrents.delay=2*ms
######################################################################
model_string=build_model_string(CA3_str,strength)
C3_c3 = Synapses(CA3, CA3,model=generalized_variable_model,pre=model_string)
C3_c3.connect(True,p=con_list[5])
C3_c3.delay=1*ms                              
############################################################################
CA1=NeuronGroup(200, eqs, threshold='v>-50*mV', reset='v=-60*mV')
CA1.v = -60*mV
model_string=build_model_string(CA1_str,strength)
SC = Synapses(CA3, CA1, model=generalized_variable_model,pre=model_string)
SC.connect(True, p=con_list[6])
###########################################################################
C1_ah=Synapses(CA1,CA1,pre=r'gi-=4*mV')
for i in range(200):
    C1_ah.connect(i,i)

'''clear up any open figures'''
cla()   # Clear axis
clf()   # Clear figure
close() # Close a figure window
'''inform user that simulation is under way'''
print 'initialized'
'''set up monitoring for each region'''
Mp = PopulationRateMonitor(Piriform)
Me = PopulationRateMonitor(Entorhinal)
Md = PopulationRateMonitor(Dentate)
Mc = PopulationRateMonitor(CA3)
Mc1 = PopulationRateMonitor(CA1)

'''example state-monitors showing to monitor values of individual cells and 
values of synapse groups see Brian documentation for more details about how to use these'''
cell_monitor= StateMonitor(CA3,'v',record=CA3)  
synapse_monitor= StateMonitor(C3_c3,'scaler',record=C3_c3[0,:])  
'''run the model'''
run(20*second)

'''average bins into firing, normalize and plot'''
rate_sump=[]
rate_sume=[]
rate_sumd=[]
rate_sumc=[]
rate_sumc1=[]
for zed in range(len(Mp.rate)/2000):
    rate_sump.append(np.sum(Mp.rate[((zed-1)*2000+5):(zed*2000+5)]/Hz)) 
    rate_sume.append(np.sum(Me.rate[((zed-1)*2000+5):(zed*2000+5)]/Hz)) 
    rate_sumd.append(np.sum(Md.rate[((zed-1)*2000+5):(zed*2000+5)]/Hz)) 
    rate_sumc.append(np.sum(Mc.rate[((zed-1)*2000+5):(zed*2000+5)]/Hz)) 
    rate_sumc1.append(np.sum(Mc1.rate[((zed-1)*2000+5):(zed*2000+5)]/Hz))
rate_sump=rate_sump/np.max(rate_sump)
rate_sume=rate_sume/np.max(rate_sume)
rate_sumd=rate_sumd/np.max(rate_sumd)
rate_sumc=rate_sumc/np.max(rate_sumc)
rate_sumc1=rate_sumc1/np.max(rate_sumc1)
plot(Mp.t[::2000]/ms, rate_sumc1,'r',Mp.t[::2000]/ms,rate_sumc,'g',Mp.t[::2000]/ms,rate_sump,'b',Mp.t[::2000]/ms,rate_sumd,'k',Mp.t[::2000]/ms,rate_sume,'c')
show()