Simulation Setup when using a List of NetCons, VecStims and Exp2Syn for multiple spike train stimulation

When Python is the interpreter, what is a good
design for the interface to the basic NEURON
concepts.

Moderator: hines

Post Reply
Liz_Stan_88
Posts: 7
Joined: Thu Jul 04, 2019 4:58 am

Simulation Setup when using a List of NetCons, VecStims and Exp2Syn for multiple spike train stimulation

Post by Liz_Stan_88 »

Good Day Mr Carnevale and Mr Hines,

I was wondering if you could assist my understanding when simulating a neuron with multiple synaptic inputs along its dendrite.
Let me describe my setup and my code to give the bigger picture, followed by the questions I have.

I am using Zilany's Model 2014 to generate a number of spike trains (i.e. you specify the number of ANFs = ANFx) for a number of center frequencies (i.e. CFx). Each generated spike train consists of an array of time instances, but it must be noted that these arrays are not of equal length/have the same number of elements per array.

The output of Zilany's Model is saved as such:

Code: Select all

		     st1,st2,st3,....,stn
		________________________________
0	C	| ANF1: [t1, t2, t3, ... , tn] |
1	F	| ANF2: [t1, t2, t3, ... , tn] |
2	:	| ANF3: [t1, t2, t3, ... , tn] |
3	1	| ANF4: [t1, t2, t3, ... , tn] |
4		| ANFx: [t1, t2, t3, ... , tn] |
		________________________________
:	C	| ANF1: [t1, t2, t3, ... , tn] |
:	F	| ANF2: [t1, t2, t3, ... , tn] |
:	:	| ANF3: [t1, t2, t3, ... , tn] |
:	2	| ANF4: [t1, t2, t3, ... , tn] |
:		| ANFx: [t1, t2, t3, ... , tn] |
		________________________________	
:	C	| ANF1: [t1, t2, t3, ... , tn] |
:	F	| ANF2: [t1, t2, t3, ... , tn] |
:	:	| ANF3: [t1, t2, t3, ... , tn] |
:	x	| ANF4: [t1, t2, t3, ... , tn] |
n		| ANFx: [t1, t2, t3, ... , tn] |
		________________________________

The reason for the explanation of the input to my NEURON model is to understand the setup of the neuron, more specifically the dendrite. The "nseg" value for the dendrite of my model is dependent on the number of centre frequencies (CF) that I choose in Zilany. The reason being is, I want each segment of the dendrite to be responsible for a single CF (i.e. I want to keep the tonotopic arrangement of the termination of the spike trains along the dendrite).

The NEURON model I created is that of a simple neuron (i.e. one soma, one dendrite, one axon) for testing purposes before I implement my own biophysically defined neuron. The code below shows that setup.

Code: Select all

"""Single Simple Neuron:"""
#Create Cell Sections:
soma = h.Section(name = 'soma')		#This is the soma.
dend1 = h.Section(name = 'dend1')		#This is the dendrite.
axon = h.Section(name = 'axon')			#This is the axon

#Define Cell Topology:
dend1.connect(soma(0), 0)			#The dendrite end 0 connects to the soma end 0.
axon.connect(soma(1),)				#The axon end (default 0) connects to the soma end 1.

#Define Geometry of Sections of the Cell:
"""Dendrite1"""
dend1.diam = 3   					#Diameter [um]
dend1.L = 250    					#Length [um]
dend1.nseg = len(CF_Range)  			#Number of Segments equals the number of CFs used in Zilany. NB must always be ODD!
"""Soma"""
soma.diam = 25  					#Diameter [um]
soma.L = 25     						#Length [um]
soma.nseg = 1   					#Number of Segments
"""Axon"""
axon.diam = 3   					#Diameter [um]
axon.L = 32     						#Length [um]
axon.nseg = 1   						#Number of Segments
h.define_shape()					

#Define the Biophysics of each Section of the Cell:
for sec in h.allsec():
    sec.Ra = 100      					# Axial resistance in Ohm * cm
    sec.cm = 0.9      					# Membrane capacitance in micro Farads / cm^2
"""Insert active Hodgkin-Huxley current in the soma:"""
soma.insert('hh')					#Inserts HH biophysics in soma.
for seg in soma:					#For all defined segments of the soma insert the following:
    seg.hh.gnabar = 0.12  				# Sodium conductance in S/cm2
    seg.hh.gkbar = 0.036  				# Potassium conductance in S/cm2
    seg.hh.gl = 0.0003    				# Leak conductance in S/cm2
    seg.hh.el = -54.3     				# Reversal potential in mV
"""Insert active Hodgkin-Huxley current in the dendrite:"""
dend1.insert('hh')					#Inserts HH biophysics in dendrite
for seg in dend1:					#For all defined segments of the dendrite insert the following:
    seg.hh.gnabar = 0.12  				# Sodium conductance in S/cm2
    seg.hh.gkbar = 0.036  				# Potassium conductance in S/cm2
    seg.hh.gl = 0.0003    				# Leak conductance in S/cm2
    seg.hh.el = -54.3     				# Reversal potential in mV
    # seg.pas.g = 0.001  				# Passive conductance in S/cm2
    # seg.pas.e = -65    				# Leak reversal potential mV
"""Insert passive current in the axon:"""
axon.insert('pas')					#Inserts pas biophysics in axon
for seg in axon:						# For all defined segments of the axon insert the following:
    seg.pas.g = 0.001  				# Passive conductance in S/cm2
    seg.pas.e = -65    					# Leak reversal potential mV

h.topology()						# Shows topology of created neuron
ps = h.PlotShape()                                     # False tells h.PlotShape not to use NEURON's gui. Shows/draws the generated Octopus Cell.
ps.show(0)
Seeing that the number of CFs can change per user input, I decided to calculate a Placement array (i.e. where each element in the array corresponds to the place along the dendrite at which that ANF for a specified CF would terminate). I kept in mind that "nseg" should be kept as an odd number, and thus ensured that the number of CFs specified should ALWAYS be odd.

My calculation of the placement array is as follows:
  • nseg = number of CFs chosen (ALWAYS ODD)
  • L = Length of each segment = 1/nseg
  • HLen = Half the length of each segment = 0.5 * (1/nseg)
  • Cfperseg = Number of segments per CF = nseg/number of CFs chosen
  • CFs chosen = [CF1, CF2, ..., CFn]
  • Placement = HLen + (x - 1)*(L), where x is the index in the list called CFs chosen
The NEURON simulation setup I use, evolved from my understanding thus far when using Netcon(), Exp2Syn() and VecStim() as list objects.

First step that I took was to restructure the output from the Zilany Model by taking into consideration that all spike trains related to a specific center frequency, which in turn related to a specific place/segment along the dendrite. Secondly, I also considered that I may use more than one dendrite, thus I may have to split the number of spike trains per center frequency evenly among the number of dendrites I specify/create. At the moment the model uses only one dendrite, but this is a future consideration. Therefore I ensured that the number of ANFs I specify in the Zilany Model is ALWAYS a multiple of the number of dendrites I create in NEURON.

Based on the above considerations I restructured the Zilany spiketrains such that I group all time elements per array index together, if the index does not exist in the original Zilany spike train output (i.e. as mentioned above each spike train does not contain the same number of elements), I append a value of 0. Therefore my spike train data I use for my NEURON model looks like this (i.e. a python List of Lists of Lists):

Code: Select all

Place on Dendrite:		|			1st nsegment on each dendrite = CF1			||				2nd nsegment on each dendrite = CF2		||				nth nsegment on each dendrite = CFn		|	
Dendrite number:		|	|	D1		|	...	|	nth D		|	||	|	D1		|	...	|	nth D		|	||	|	D1		|	...	|	nth D		|	|
0 	Simrun 1 = st1	[ 	[	[t1, t1, ...... ,nth t1	] 	...	[t1, t1, ..... , nth t1	]	][	[t1, t1, ...... ,nth t1	] 	...	[t1, t1, ..... , nth t1	]	][	[t1, t1, ...... ,nth t1	] 	...	[t1, t1, ..... , nth t1	]	]	]
1 	Simrun 2 = st2	[ 	[	[t2, t2, ...... ,nth t2	] 	...	[t2, t2, ..... , nth t2	]	][	[t2, t2, ...... ,nth t2	] 	...	[t2, t2, ..... , nth t2	]	][	[t2, t2, ...... ,nth t2	] 	...	[t2, t2, ..... , nth t2	]	]	]
: 	   :		: 	:			:				:					:					:					:					:									
n 	Simrun n = stn	[ 	[	[tn, tn, ...... ,nth tn	] 	...	[tn, tn, ..... , nth tn	]	][	[tn, tn, ...... ,nth tn	] 	...	[tn, tn, ..... , nth tn	]	][	[tn, tn, ...... ,nth tn	] 	...	[tn, tn, ..... , nth tn	]	]	]

This follows the logic I used for the setup of the NEURON simulation as follows:
  • For each segment along the dendrite to which I intend to input a synaptic input, I create an h.Exp2Syn() and append that to an h.List(). Thus each segment along the dendrite would have its own unique h.Exp2Syn() point process.
  • I also create one h.VecStim() for each time instance per segment along each dendrite. I also append each h.VecStim() to a h.List()
  • For each h.VecStim() created I create an h.Vector() to which the specified spike time per dendrite per segment of dendrite is inserted and that h.VecStim() can play.
  • I then create a h.NetCon() object for all the h.VecStim() using the h.Exp2Syn() for that specific segment of the dendrite.
  • I then set up the recording vectors for the tvec and voltage membrane of the soma.
  • Finally, for each simrun I call h.run() and then I reinitialise repeating steps 1 to 4 before h.run() again.
My code for this setup described above is given below.

Code: Select all

"""NEURON Simulation Run Setup:"""
#++++++++++++++++++++++++++++++++++++++
#Simulation Environmental Constants:
#++++++++++++++++++++++++++++++++++++++
tau1 = 0.07                                                     	#DECAY TIME IN ms.
tau2 = 0.34                                                     	#RISE TIME IN ms.
e = 0.0                                                      	#REVERSAL POTENTIAL IN mV.
weight = 2							#NetCon Weight [nS].
delay = 0								#NetCon Delay = 0.
threshold = 0							#NetCon Threshold = 0.
Temp = 34                                                       	#TEMPERATURE OF ENVIRONMENT.
T_Step = 0.005                                                  #TIME STEP.
v0 = -63                                                        	#INITIAL VOLTAGE MEMBRANE IN mV.
Tstop = 10                                                     	#TIME OF STIMULATION IN s.
V_init = -63                                                    	#Initial voltage in mV of membrane.
Vm_soma = []							#List to append Vm Soma results to per simrun.
tvec = []								#List to append time vec results to per simrun.
Syn_I = []								#List to append Synaptic Current results to per simrun.
ALL_synI = []							#List to append ALL synaptic current results to.
for simrun in simrunList:					#Per simrun/ number of elements in simrunList (A list of lists) execute the following:
    h.celsius = Temp                                      	#Sets the NEURON Environment temperature
    h.dt = T_Step                                                #Sets the NEURON stimulation time step
    h.v_init = V_init                                            #Sets the NEURON initial voltage
    h.finitialize(v0)                                    		#Initializes the NEURON initial voltage
    exp2synList = h.List()					#Creates a h.list to append all exp2syn created.
    vecstimList = h.List()					#Creates a h.list to append all vecstim created.
    netconList = h.List()					#Creates a h.list to append all netcons created.
    synIList = h.List()						#Creates a h.list to append all synI created.

    for placealongdend in simrun:			#For each segment along the dendrite per simrun execute the following:
        dend = 1
        for dendnum in placealongdend:		#For each dendrite at the specified segment per simrun execute the following:
            if dend == 1:							#If the dend num = dendrite 1 execute the following.
                dendtocall = dend1(Placement[place])		#Create variable dendtocall which specifies the neuron section and segment
            elif dend == 2:						#else If the dend num = dendrite 2 execute the following.
                dendtocall = dend2(Placement[place])		#Create variable dendtocall which specifies the neuron section and segment
            elif dend == 3:						#else If the dend num = dendrite 3 execute the following.
                dendtocall = dend3(Placement[place])		#Create variable dendtocall which specifies the neuron section and segment
            elif dend == 4:						#else If the dend num = dendrite 4 execute the following.
                dendtocall = dend4(Placement[place])		#Create variable dendtocall which specifies the neuron section and segment
            else:								#else execute the following.
                print("Number of dendrites != num of dends specified! CHECK CODE")
                pass
                
            """setup Exp2Syn"""					#The EXP2SYN setup
            synapse = h.Exp2Syn(dendtocall)			#Create h.Exp2Syn() point process for specific dendrite at specific dendrite segment.
            exp2synList.append(synapse)				#Append synapse to exp2synList.
            exp2synList[-1].tau1 = tau1				#Set decay time.
            exp2synList[-1].tau2 = tau2				#Set rise time.
            exp2synList[-1].e = e					#Set reversal potential.

            """Setup Vecstim"""						#The VECSTIM setup
            for time in dendnum:					#For each time element in list execute the following:
                vecstim = h.VecStim()					#Create h.VecStim() for time element.
                vecstimList.append(vecstim)				#Append vecstim to vecstimList.
                spiketime = h.Vector([time])				#Create a h.Vector() and insert the time element in the vector.
                vecstimList[-1].play(spiketime)				#Play the spiketime.

                """Setup Netcon:"""					#The NETCON setup	
                net = h.NetCon(vecstimList[-1], exp2synList[-1])	#Create h.NetCon() using same exp2syn for all time elements in List per dendrite per segment of dendrite.
                netconList.append(net)						#Append net to netconList.
                netconList[-1].weight[0] = weight				#Set Netcon weight.
                netconList[-1].delay = delay					#Set Netcon delay.
                netconList[-1].threshold = threshold				#Set Netcon threshold.

            dend += 1									#Increment dendrite number.

    """Run Simulation"""
    somaR = h.Vector()						#Create a vector for soma results.					
    somaR.record(soma(0.5)._ref_v)				#Record the vm of the soma.
    timeR = h.Vector()						#Create a vector for time results.
    timeR.record(h._ref_t)						#Record the simulation time.
    for es in exp2synList:						#For each exp2syn created execute the followiing:
        synI = h.Vector()							#Create a vector for the synaptic current results
        synIList.append(synI)							#Append the synaptic current vectors to synIList
        synIList[-1].record(es._ref_i)					#Record the synaptic current.

    h.tstop = Tstop							#Set the simulation duration.
    h.run()									#Run NEURON simulation.

    Vm_soma.append(np.array(somaR))			#Append Soma results to Vm_soma List
    tvec.append(np.array(timeR))				#Append time results to tvec List
    for sI in synIList:							#For all synapse results in synIList execute the following:
        Syn_I.append(np.array(sI))					#Append all synaptic current results to Syn_I List
    ALL_synI.append(Syn_I)						#Append Syn_I to ALL_synI List.

print("simulation Done")					#print to cmd simulation is done.
My questions or concerns, arise as I am the only member in my research group that has ever used NEURON software, so I have definitely been learning by error and learning from the Forum. I guess all good skills to learn being in research and persuing the career of researcher, but my problem is I do not have a fellow or well versed in NEURON colleague to bounce ideas or concerns off of. Thus, any help or guidance will be greatly appreciated. My questions/concerns are as follows:
  • Seeing that I am only interested in the interspike intervals of the soma response, should I be inserting a h.NetCon() to record the tvec at the soma?
  • Should I be re-initialising the parameters per Simrun, as seen in the above code? I ask, because in my mind by reinitializing I make all the parameters in the simulation = 0 or their equivalent initial values, which would basically mean a new cell per run. Am I correct? As in real life all these spike times would be arriving at the same cell, however consecutive action potentials will produce if threshold at the soma is exceeded and this is also dependent on how fast/slow the soma's ability is to polarize and depolarize based on its ionic channels.
  • Is there not a way to cut out the restructuring of the Zilany Model Output, i.e. h.Vector([each spike train per CF]) or am I correct in my setup and use of h.VecStim and h.Vector, seeing that timing characteristics are important? Trying to improve computational efficiency I guess.
  • Is appending 0 when restructuring the Zilany model as I mentioned above not influencing my data as an Exp2Syn() would occur at time t = 0, when in actuality no synaptic input or arrival of spike time on the specific segment of the dendrite should occur? Thus, is there not a NEURON shorthand telling NetCon/VecStim not to produce an event if the spike time = 0? The only way I have thought about this is simply implementing an "IF" "ELSE" statement to basically skip over the creation of a h. VecStim() if t = 0.
  • Does NetCon.weight not need to vary along the dendrite length? I guess I am asking cause I understand that the weight variable measured in nS is the weight or for better words the strength of the Exp2Syn point process. Thus the further away from the soma the more weight should be applied to produce an action potential of same strength to an action potential produce by a synaptic input closer to the soma. Or is this totally not necessary.
  • Furthermore, eventhough I have read up on h.NetCon(), I do not fully understand Netcon.threshold and how this influences the synaptic input. A little bit of enlightenment would be appreciated?
  • Lastly, I think my biggest concern is, IS the synaptic instance along the dendrite being applied as I expect it should (seeing I'm newish in the use of NEURON and definitely new to coding using h.List())? Namely, have I forgotten/ not forgotten to clear a pointer or perhaps forgotten to consider something that will potentially affect the production of action potentials at my Soma? I do believe, that my output of my created model to be correct based on my results when I do plot my synaptic currents and voltage membrane potential of soma. Furthermore, I have checked that for each NetCon() object I am using the correct segment, source and target, as well as checking that each Exp2Syn() point process that I applied has the specified tau1, tau2, e and what the current it would be.
In conclusion, I do apologise for the extremely long query, however, I have learnt that being detailed tends to help the person trying to assist by giving them enough information to see where the problem lies. Again any assistance or guidance will be much appreciated.

Regards
ramcdougal
Posts: 224
Joined: Fri Nov 28, 2008 3:38 pm
Location: Yale School of Public Health

Re: Simulation Setup when using a List of NetCons, VecStims and Exp2Syn for multiple spike train stimulation

Post by ramcdougal »

I haven't wrapped my head around all of this, but here's a few points to help:

I'm not sure I understand your placement array, but there's a comment "L = Length of each segment = 1/nseg"... the length of a segment is the length of the section (which presumably isn't 1 micron) divided by that section's nseg. So for example: dend_seg_length = dend.L / dend.nseg

There's no advantage to using h.List here; a regular Python list would work just fine. i.e. exp2synList = [] and so forth

A NetCon's threshold has to do with when it decides the presynaptic cell has fired (the membrane potential crosses the threshold). It matters for recording output spike times (in your case from the soma). It has no effect with NetStims or VecStims that are told explicitly when the stimulus is to happen (and thus have no presynaptic cell with no membrane potential to monitor). If it's true that all you care about are ISIs, then you could just have a NetCon at the soma that records spike times and then you don't need Vm_soma or tvec. (But: since you have them, you might as well test the NetCon spike times by making sure they line up with the time series recording. Here's an example: https://colab.research.google.com/drive ... sp=sharing

Note that `h.run()` does two conceptual things: it reinitializes your simulation (not necessarily setting to 0 but applying whatever initialization rules have been specified) and then runs the simulation to h.tstop. In particular, this means that your explicit call to h.finitialize(v0) earlier has no effect. There are some edge cases (with legacy code) where the following won't work, but I recommend for most cases including yours, to initialize explicitly and to continuerun to a specific time (this assumes initialization has already been done); e.g.

Code: Select all

from neuron import h
from neuron.units import mV, ms

# setup model here

h.finitialize(-65 * mV)
h.continuerun(20 * ms)
Continuing from that, as h.run() begins by initializing the simulation, you see that this in particular does not break the structure of the model. That is, you could have created e.g. the exp2synList outside of the for loop. It's the states that change in initialization.

Depending on how long this series of runs takes, you may want to consider using a database (Python comes with sqlite3) and dumping the data there as it is generated so that you can look at it and run some tests while the simulations are continuing.

Finally a Python comment: you never need to manually advance an index like you're doing with dend. If you really need the index, you can get it with enumerate, but here you don't need the index; you just need to select the dendrite. So in other words, instead of incrementing dend and using a big if statement, you could just:

Code: Select all

for my_dend, dendnum in zip([dend1, dend2, dend3, dend4], placealongdend):
    dendtocall = my_dend(Placement[place])
Liz_Stan_88
Posts: 7
Joined: Thu Jul 04, 2019 4:58 am

Re: Simulation Setup when using a List of NetCons, VecStims and Exp2Syn for multiple spike train stimulation

Post by Liz_Stan_88 »

Thank you so much ramcdougal for your quick response. I do understand I put alot of information in my query and it definitely is a lot to get one's head wrapped around.
I'm not sure I understand your placement array, but there's a comment "L = Length of each segment = 1/nseg"... the length of a segment is the length of the section (which presumably isn't 1 micron) divided by that section's nseg. So for example: dend_seg_length = dend.L / dend.nseg
With regards to my placement calculation, I figured out today I do not need to calculate one, to be able to point to the correct segment for the dendrite. But the reason I did was based on my understanding of how NEURON deals with the cable equation, so although you specify a specific length for the dendrite in this case 250 um, when you apply an nseg > 1 it equivalently as you state above divides that length up into the number of segments specified by nseg. However this length you have calculated [in um] is not the segment value or point to which NEURON would point to. Because the entire section is converted to a value between 0--(nsegs)--1. Hence why the length of my segment calculation is 1/nseg. For future users having the same issue when trying to point to a specific segment in your section the following code does what my calculation does in two lines of code, where sec is the defined section you want to for loop over in my case it would be dend1 and x is the value of the segment:

Code: Select all

for seg in sec:
    print('{}'.format(seg.x))

There's no advantage to using h.List here; a regular Python list would work just fine. i.e. exp2synList = [] and so forth
I thank you for this clarification I was originally using a regular Python list, so I will revert back to its use.
A NetCon's threshold has to do with when it decides the presynaptic cell has fired (the membrane potential crosses the threshold). It matters for recording output spike times (in your case from the soma). It has no effect with NetStims or VecStims that are told explicitly when the stimulus is to happen (and thus have no presynaptic cell with no membrane potential to monitor). If it's true that all you care about are ISIs, then you could just have a NetCon at the soma that records spike times and then you don't need Vm_soma or tvec. (But: since you have them, you might as well test the NetCon spike times by making sure they line up with the time series recording. Here's an example: https://colab.research.google.com/drive ... sp=sharing
Note that `h.run()` does two conceptual things: it reinitializes your simulation (not necessarily setting to 0 but applying whatever initialization rules have been specified) and then runs the simulation to h.tstop. In particular, this means that your explicit call to h.finitialize(v0) earlier has no effect. There are some edge cases (with legacy code) where the following won't work, but I recommend for most cases including yours, to initialize explicitly and to continuerun to a specific time (this assumes initialization has already been done); e.g.
Thank you for this clarification, it now makes more sense and I think I will add a h.Netcon() at my soma to confirm that the spike times line up as you suggested. Furthermore, I will adjust my code using the h.continuerun(), as this makes more logical sense to me as well, I just was really unsure how to implement it within the massive for loop structure I had going on.

I will also be using your python (sqlite3) suggestion as although the run takes approximately 45min in total to complete, I could use that 45 min more productively than sitting around waiting for everything to finish.

Finally your last suggest, helps me tremendously, I was actually searching for a solution to the "if", "else" statement for getting rid of manually advancing my index for dend. So this was definitely a win.
ted
Site Admin
Posts: 5885
Joined: Wed May 18, 2005 4:50 pm
Location: Yale University School of Medicine
Contact:

Re: Simulation Setup when using a List of NetCons, VecStims and Exp2Syn for multiple spike train stimulation

Post by ted »

Comments:

First, your code is still very much in the development and debugging stage. At this point there is no particular advantage to running long simulations, or even tens of short simulations. Start with a single simulation that lasts only long enough for one spike train--or even a quarter of a spike train. Only when the model's response to that makes sense will it be time to increase simulation duration and/or the number of runs.

But before even that, your primary concern at this stage is to verify that there is a close match between your conceptual model and what is in the computer. In this case, "conceptual model" includes anatomical and biophysical properties of the model cell, and the nature of the stimuli (their temporal structure and where they are located in the cell). And, since this implementation draws on a different model by somebody else, are you concerned about the match or mismatch between that previous model (which was very much a signal processing model--no neurons or synapses in it) and your current conceptual model?

Second, just recording "spiking" at the soma is not enough. During development you should also examine the time course of your model cell's membrane potential. How else can you know whether your model cell is actually "spiking" or whether it is simply being forced above NetCon.threshold by synaptic inputs that are nonphysiologically large and prolonged? Furthermore, the Hodgkin-Huxley spike mechanism is only capable of spiking over a rather limited range of frequencies, and it has a wide spike that is followed by even wider absolute and relative refractory periods. That imposes limitations on the kinds of spike trains that your model cell will be able to respond to without dropping spikes.

With regard to the placement of synapses, it isn't clear to me that your strategy has anything to do with the Zilany model, but that could be a good thing if all you want from their model is the temporal structure of inputs. So let's assume it's a good thing--you just want a tonotopic arrangement of synapses along one or more dendrites. I would suggest that instead of thinking "one synapse per segment" you should think either "one synapse per anatomical distance dx along a section" or maybe better "one synapse per normalized distance "drange" (as in delta range) along a section." The latter is easily done. If you have N synapses to distribute evenly along a section, forget about the section's nseg parameter. Instead do this (in python-tinged pseudocode)

Code: Select all

for x in range(N):
  attach an ExpSyn to sectionname((x+0.5)/N)
  append the ExpSyn to a list
This takes advantage of the fact that
for x in range(N)
iterates from 0 to N-1
and
(x+0.5)/N
converts that sequence to N internal locations that are evenly spaced along the section from 0.5/N at one end to (N-0.5)/N at the other--and if N is odd, there will always be one location equal to 0.5.

Notice that nseg doesn't appear anywhere in the calculation. Instead, you can specify nseg to have any value you like (e.g. a value sufficiently large to get good spatial accuracy while avoiding unnecessary computational burden).

"Well, if nseg is too small, multiple synapses might be attached to the same segment. I really want each segment to deal with its own frequency."

It won't. In a properly discretized cable (i.e. sufficient for computational accuracy) the coupling between adjacent compartments is so tight that no segment is going to "see" only the current delilvered by its own attached synapse. Charge will redistribute to adjacent segments in a very few ms, which is much faster than the minimum interspike intervals that the HH mechanism can generate.

But you can see this for yourself if you run a few simulations with your model cell(s) and look at variables other than (apparent) "soma spike times."
Post Reply