Running secondary networks in isolation

Moderator: wwlytton

Post Reply
psipeter
Posts: 7
Joined: Tue Mar 28, 2017 5:12 pm

Running secondary networks in isolation

Post by psipeter »

I am running NEURON with Python within the context of another neural simulator, Nengo. The ordering of the build process in Nengo requires that I create all the NEURON objects (cells from hoc template, NetCons, etc.) at an early stage. Later in the build process, I temporarily create isolated populations of cells, simulate them, collect data, and delete them. This data helps me optimize high-level Nengo parameters (connection weights), but does not require any connection to, or information about, the original NEURON objects.

My problem is that, in simulating these isolated populations, I am also inadvertently simulating all the objects I created earlier in the build. Although these processes do not interfere with the isolated populations or corrupt the data, NEURON continues to simulate the cells, drastically and unnecessarily increasing runtime.

My question: is it possible to either (a) create an isolated network populated with test cells and simulate only this network, even after cells have been created elsewhere in the same thread; or (b) create all the NEURON objects, then set a flag on some of them that tells the simulator to disregard them while other objects are being run?

Thanks,
- Peter
hines
Site Admin
Posts: 1682
Joined: Wed May 18, 2005 3:32 pm

Re: Running secondary networks in isolation

Post by hines »

It is not possible to not simulate created cells in NEURON. I can imagine this is a problem if the expense of setting up a network takes a long
time. Otherwise, the best way to proceed is to destroy the previous network and setup only a network needed for the next phase. This can
be most easily done by encapsulating the network setup in a class with an appropriately parameterized constructor. Then to destroy a network
it suffices to make the instantiated class reference count go to 0. By the way, when destroying a network, it is best to first clear the gids (if
used) with pc.gid_clear(), then destroy the list of NetCon, then destroy the list of cells.

Does your problem reasonably fit into the style of successively building and tearing down networks?
psipeter
Posts: 7
Joined: Tue Mar 28, 2017 5:12 pm

Re: Running secondary networks in isolation

Post by psipeter »

Thanks for the rapid response. I had previously considered the build/teardown approach, and in fact I use it for a preliminary optimization elsewhere in my code. However, the ordering of the Nengo build process makes this quite problematic in the current context.

In Nengo, each ensemble is populated early in the build (in succession), then later each connection between ensembles is constructed (in succession) using a function called build_connection(). This function is the only place in the build process where I have access to the parameters I need to run my optimization, and hence this is where I perform the temporary creation/teardown of a stand-alone NEURON network. However, build_connection() is also the logical place to construct the NetCon objects that will be used in the simulation proper. If I create these objects after optimizing one connection between ensembles, then build_connection() gets called again to form a second connection, the second optimization will unnecessarily simulate the cells/NetCons created during the first call.

The only solution I can see is to wait to actually create the NEURON cells/NetCons until after all the builder has finished all calls to build_connection(). It may be possible to hack such a method into the build order, but it certainly would not respect the overall Nengo build flow, and would likely make my NEURON interface incompatible with the core Nengo release. I’ll consult some of the other developers.
hines
Site Admin
Posts: 1682
Joined: Wed May 18, 2005 3:32 pm

Re: Running secondary networks in isolation

Post by hines »

Sorry. I'm not quite seeing the problem yet. I'm hearing that when build_connection is called, that you need to build something that you can
optimize and then throw it away. Can you not also save the relevant info of each build_connection call (and return value) in set/dict/list for later simulation proper? It may be that is the workflow of the "only solution" you mention. I gather from the "However...unnecessarily simulate... created
during the first call" sentence that the abstract info you need to save might be unreasonably large?...or take too much time to keep building the
same portions of the network too many times for reasonable setup performance? Would you prefer to be able to mark Sections (including all their
synapses and mechanisms) or
ARTIFICIAL_CELL instances as "Do not Simulate"?

Just curious. is parallelization important in the simulations?
psipeter
Posts: 7
Joined: Tue Mar 28, 2017 5:12 pm

Re: Running secondary networks in isolation

Post by psipeter »

You’re correct that the solution I’m proposing is to save the relevant optimized info (connection weights) for later use in the simulation proper. The problem with this solution isn’t that this information is too large, or that creation/destruction of portions of the network is too expensive. It’s that, by the time build_connection() is called, Nengo has already created all the NEURON sections for the simulation proper. It is difficult or impossible to delete all these previously-built cells (before the first build_connection() call) then recreate them later with the optimized info (after the last build_connection() call) because of the limited control I have over Nengo’s build order. Does that make more sense?

From my perspective (working within the Nengo environment), a less-intrusive solution would be to mark Sections with “Do not Simulate” during the build process, then switch them to “Do Simulate” after all the optimization is complete. And no, parallelization is not important.
hines
Site Admin
Posts: 1682
Joined: Wed May 18, 2005 3:32 pm

Re: Running secondary networks in isolation

Post by hines »

I was completely unaware that Nengo can create a NEURON Section. I don't know very much about Nengo but thought it had its own native
implementation of models. Can you point me to a document that explains more of the Nengo workflow related to NEURON.
As I think about the possibility of removing a Section from the list of Sections to be simulated, I keep remembering implementation details that
make it more and more remote as a possibility (mostly due to at least two implementation layers between description and the structures involved
in the vectorized number crunching. Nevertheless, it is interesting for me to consider what it would mean for a Section to exist in a hoc_List other
than the internal
extern hoc_List* section_list; /* Where the Sections live */
I don't know if that would support a <donotsimulate> flag or if everything would collapse in a pile of rubble. It might just work if a Section had
no path to it from a simulated Section i.e every section of a tree cell would have to be moved to the "do not simulate" list.
psipeter
Posts: 7
Joined: Tue Mar 28, 2017 5:12 pm

Re: Running secondary networks in isolation

Post by psipeter »

Haha, well you’re probably unaware of the nengo-NEURON interface because the only one working on this project is...me :) Let me actually introduce myself - I’m Peter Duggins (psipeter@gmail.com), a 2nd year PhD student in Chris Eliasmith’s lab. As part of my masters thesis (also with Chris), I developed some code that allows Nengo to build and simulate small populations of NEURON objects together with traditional Nengo neuron_types (typically spiking LIF). Although I managed to get these networks working and have applied the associated optimization methods to build models of working memory (in part to investigate the effects of biophysical perturbation, corresponding to the injection of pharmacological agents, on short-term memory), the interface is far from perfect, and continues to be under development. As such, there isn’t a single document I can point you to that would be helpful in this regard. I’d be more than happy to share with you the work-in-progress code, go through some of the workflow and see if we can find non-invasive solutions, but I’m reluctant to post the github link on a forum until I’ve got a stable release and a nice written set of instructions.

Having worked with Nengo and Neuron for a few years now, the biggest incompatibilities I’ve encountered between the two platforms result from the issue we’re discussing here: not being able to simulate a subset of NEURON objects. Although this isn’t a problem during the main simulations, there are several points during the build process when Nengo uses various principles from the Neural Engineering Framework to parameterize the network. Without going into a full explanation of the NEF, Nengo tries to assign each neuron a “gain” (which affects the scaling of spiking inputs to the cell) and a “bias” (which represents the background input to the cell) in a way that creates a diverse set of “tuning curves” (firing rate vs. vector-valued state-space input) among a neural population. Normally, when Nengo solves for gains and biases, it repeatedly runs a neuron_type’s step_math() method (as you said, native implementation of the state update equations) with various current inputs. The idea is to approximate the neurons’ responses without actually simulating the whole network. While my version of this optimizations is a bit more advanced, it still relies on the assumption that I can isolate a group of neurons, simulate them to determine their behavior, add them to the network proper, then move on to optimize the next population. While I can do all these things, networks which include more than one population of/connections to NEURON objects will slow down with each successive optimization.

I can appreciate how difficult it could be to prevent Sections from being simulated, given that NEURON is designed around efficiently performing simultaneous vectorized computations. Not trying to make your life miserable just to add one feature that breaks everything else. But if you think it’s at all possible, it would really facilitate the merger of these two simulations platforms. In my biased opinion, NEURON and Nengo respectively represent the best low- and high-level neural modelling languages out there, and combining them opens up many new avenues in unifying biophysical and functional modelling for neuroscientists.
hines
Site Admin
Posts: 1682
Joined: Wed May 18, 2005 3:32 pm

Re: Running secondary networks in isolation

Post by hines »

You've convinced me that removing cells from a simulation without
destroying them is a useful feature. I've concluded that the basic
functionality will have a straightforward implementation. More
rarely used combinations of features that will break the system
can raise errors until the implementation is extended to encompass those
features. Examples of such features are use of
cvode.cache_efficient(1) and the use of gap junctions or any pointers
in one cell to a variable of another. I suppose gap junctions will work
in a formal sense in that the non-simulated target will represent a
constant voltage but it is likely best to set that gap junction conductance
to 0 until the target is reactivated. Similarly, NetCon events, unless
NetCon.active(0), will be processed by the target synapse NET_RECEIVE block.
That is not really a problem as it merely results in some unused state
variables changing discontinuously. There are likely a few things I have
not considered, but I don't see a problem with variable time step methods
or parallel simulations. Even threads should be ok modulo the
cvode.cache_efficient(1) wildcard which might require a rewrite to bring
it up to date with the subsequent introduction of NMODL pointer semantics
information. I believe almost all spike coupled network simulation would
be encompassed by the basic functionality implementation.
Validation can consist in results for the simulated subset being the same
as if the whole system is simulated.

I will pursue this in greater depth and get back to you reasonably soon
with a "nosimulate" branch in the github.com/neuronsimulator/nrn
repository. Or an apology that it was sadly not feasible. In the
former case, you can then compare time and results in your Nengo context.
It should suffice at the user level to make available in Python a
nrnSection.simulate(bool) # can use 0 or 1 as well
which, with an arg of False, moves the currently accessed
section and all the sections in its tree (cell) to the list of sections
not to be simulated. With no argument, it will return True or False.
Note that every variable associated with a "not simulated"
cell can be assigned and evaluated at the interpreter level.
The Python idiom would have a form like
for cell in celllist:
cell.soma.simulate(0)
but let me know if you need a hoc interpreter version in which case I'd
make a cvode.do_not_simulate(1) method which in python would be used as
h.cvode.do_not_simulate(1, sec=cell.soma)
psipeter
Posts: 7
Joined: Tue Mar 28, 2017 5:12 pm

Re: Running secondary networks in isolation

Post by psipeter »

Fantastic, thanks so much for taking a look at this! The python idiom you sketched would be perfect for my purposes. Keep me posted on updates, and I'll do the same for the Nengo implementation if/when the NoSimulate feature becomes active.
hines
Site Admin
Posts: 1682
Joined: Wed May 18, 2005 3:32 pm

Re: Running secondary networks in isolation

Post by hines »

Ok. The do-not-simulate branch is ready to try. the log message is
commit 0ebc0f8ff10ddd68b8cd3fe12cfb5e78bd9dd6a5 (HEAD -> do-not-simulate, origin/do-not-simulate)
Author: Michael Hines <michael.hines@yale.edu>
Date: Thu Jan 24 13:38:54 2019 -0500

nrn.Section.simulate(False) causes the entire tree containing the section
to be removed from the internal list of sections to be simulated.
nrn.Section.simulate() will return True if the section is being simulated.
The relevant POINT_PROCESS and density mechanisms will also not be simulated.
(although NetCon delivery will proceed normally).
Note that unsimulated sections will not appear in the "for sec in h.allsec():"
loops, h.topology(), etc. So the user should maintain references to the
unsimulated sections in order to be able to put them back into the internal
list of simulated sections via nrn.Section.simulate(True). Alternatively
one may use section.simulate(-1) to force all sections to be simulated.
Note that if there are unsimulated sections, cvode.cache_efficient() is
turned off and refuses to be turned on until the number of unsimulated sections
is 0.
I tested it minimally with a series of 5 rings each with 8 cells with each ring successively created and simulated so that a spike went around
the ring several times. Then the entire ring was unsimulated and the next ring was started. At the end, all rings were simulated. For each
simulation the spikes were recorded and appended to a file. I experienced two curiosities. First, I had to inactivate all NetCon to cells that were
unsimulated to prevent "NEURON: ExpSyn[0] :Event arrived out of order.". I guess it was possible for a NetCon spike to be in transit which then
got delivered when the simulation started running the next ring. (unsinulated cells don't get initialized). Second I had to set the voltage of the soma of each cell below its
threshold to prevent the spike dectector from generating an initial spurious spike in the case were the previous simulation happened to have a cell
above threshold when the simulation ended.
The idea is...

Code: Select all

def runs():
  types = range(0, ncell*nring)
  ncs = h.List("NetCon") 
  for i in range(nring):
    rings.append(Ring(ncell, nbranch, ncompart, i*ncell, types))
    randomize(rings[i])
    run(tstop)
    for cell in rings[i].cells:
      cell.soma.simulate(0)
      cell.soma.v = -65 #prevent generation of spike
    for nc in ncs: 
      nc.active(0) 

  rings[0].cells[0].soma.simulate(-1)
  for nc in ncs: 
    nc.active(1) 
  run(tstop)
As far as I know, this can be used with threads and mpi parallel. Only cache_efficiency is not feasible at present without a lot more work.

PS. I realized the two curiosites could be avoided by doing an h.finitialize(-65) after each run (after the spikes were saved).
psipeter
Posts: 7
Joined: Tue Mar 28, 2017 5:12 pm

Re: Running secondary networks in isolation

Post by psipeter »

Thanks for the quick patch. I'm having some trouble iterating through all the sections in a cell and setting sec.simulate(0) for each of them.

This code works fine:

Code: Select all

import neuron
import numpy as np

soma = neuron.h.Section()
dend = neuron.h.Section()
dend.connect(soma(1))
v_rec = neuron.h.Vector()
v_rec.record(soma(0.5)._ref_v)
stim = neuron.h.IClamp(soma(0.5))
stim.dur = 10  # ms
stim.amp = 10  # namp

soma.simulate(0)
dend.simulate(0)
    
neuron.h.run()

print(np.array(v_rec))
Accessing the sections with neuron.h.allsec() fails with a segmentation fault:

Code: Select all

import neuron
import numpy as np

soma = neuron.h.Section()
dend = neuron.h.Section()
dend.connect(soma(1))
v_rec = neuron.h.Vector()
v_rec.record(soma(0.5)._ref_v)
stim = neuron.h.IClamp(soma(0.5))
stim.dur = 100  # ms
stim.amp = 10  # namp

for sec in neuron.h.allsec():
    sec.simulate(0)
    
neuron.h.run(100)

print(np.array(v_rec))
Is there a different way I should be accessing the list of sections in a cell model, so that I can save each to a python list/dictionary to be re-activated later?
hines
Site Admin
Posts: 1682
Joined: Wed May 18, 2005 3:32 pm

Re: Running secondary networks in isolation

Post by hines »

You only need to call sec.simulate(0) for one section in a cell and all the sections in that cell will follow along. However, if you call
sec.simulate(0) for several or all of the sections in a cell it should not matter as any after the first will do nothing.

Code: Select all

for sec in neuron.h.allsec():
    sec.simulate(0)
It appears that the allsec iterator was not implemented for the case where sections are disappearing in the middle of the iteration. Til I fix that
a work around is

Code: Select all

for sec in [sec for sec in neuron.h.allsec()]:
    sec.simulate(0)
But to run your eniire fragment successfully you will also need either
from neuron import gui
or
neuron.h.load_file("stdrun.hoc")
Post Reply