Parallelizing sections within a single cell in python without using GUI

General issues of interest both for network and
individual cell parallelization.

Moderator: hines

Post Reply
AtulK
Posts: 1
Joined: Wed Mar 23, 2022 7:22 pm

Parallelizing sections within a single cell in python without using GUI

Post by AtulK »

As a proof-of-concept prior to applying to single cell with more than 2500 compartments, we started with parallelizing a single cell with 3-sections - soma, axon1 and axon2 to understand the process. Here is the code:

from math import pi
from neuron import h
h.load_file('stdrun.hoc')
#h.load_file("parcom.hoc")
import multiprocessing


# Simulation parameters
h.dt = 0.025 # time step (resolution) of the simulation in ms
h.v_init= -50 # initial membrane potential in mV

# Create the soma section and define the parameters
soma = h.Section(name='soma')
soma.diam = 1000; # micrometers
soma.L = 1e5/pi # micrometers
soma.cm = 9e-3 # membrane capacitance uF/cm2
soma.Ra = 10000 # ohm-cm

soma.insert('leak'); soma.insert('cas'); soma.insert('cat')
soma.insert('kA'); soma.insert('kca'); soma.insert('kdr')
soma.insert('capool')
soma.el_leak=-50; soma.eca = 120; soma.ek = -80 # mV
soma.cainf_capool = 500e-6; soma.casi = soma.cainf_capool # mM

# Create the axon section and define the parameters
axon1 = h.Section(name='axon1')
axon1.diam = 1000 # micrometers
axon1.L = 1e5/pi # micrometers
axon1.cm = 1.5e-3 # uF/cm2
axon1.Ra = 10000 # ohm-cm

axon1.insert('leak'); axon1.insert('na'); axon1.insert('kdr')
axon1.el_leak = -60; axon1.ena = 50; axon1.ek = -80 # mV

#Create second axon for parallelization
axon2 = h.Section(name='axon2')
axon2.diam = 1000 # micrometers
axon2.L = 1e5/pi # micrometers
axon2.cm = 1.5e-3 # uF/cm2
axon2.Ra = 10000 # ohm-cm

axon2.insert('leak'); axon2.insert('na'); axon2.insert('kdr')
axon2.el_leak = -60; axon2.ena = 50; axon2.ek = -80 # mV

# Connect axon to the soma
axon1.connect(soma(1),0)
axon2.connect(axon1(1),0)
h.topology()

# Current clamp
ccl = h.IClamp(soma(0.5))




h.nrnmpi_init()
pc = h.ParallelContext()
#cpu = multiprocessing.cpu_count()
#h.load_file("parcom.hoc")
#p = h.ParallelComputeTool()
def activemodel(tstop,soma_gl,gcas,gcat,gkA,gkca,soma_gkdr,taucas,axon_gl,axon_gna,axon_gkdr,amp,dur,delay):
soma.glbar_leak = soma_gl
soma.gcasbar_cas = gcas
soma.gcatbar_cat = gcat
soma.gkAbar_kA = gkA
soma.gkcabar_kca = gkca
soma.gkdrbar_kdr = soma_gkdr
soma.taucas_capool = taucas
axon1.glbar_leak = axon_gl
axon1.gnabar_na = axon_gna
axon1.gkdrbar_kdr = axon_gkdr
axon2.glbar_leak = axon_gl
axon2.gnabar_na = axon_gna
axon2.gkdrbar_kdr = axon_gkdr
ccl.amp = amp
ccl.dur = dur
ccl.delay = delay

h.tstop = tstop

pc.runworker()
#axon1.disconnect()
#axon2.disconnect()
pc.multisplit(soma(1), 1)
pc.multisplit(axon1(0), 1)
pc.multisplit(axon1(1), 2)
pc.multisplit(axon2(0), 2)
pc.multisplit()
h.run()
pc.done()
h.quit()

# default setting
tstop = 1000 # how long to run the simulation in ms
soma_gl = 0.000045e-3 # S/cm2
gcas = 0.055e-3 # S/cm2
gcat = 0.0552e-3 # S/cm2
gkA = 0.200e-3 # S/cm2
gkca = 0.500e-3 # S/cm2
soma_gkdr = 1.890e-3 # S/cm2
taucas = 303 # ms
axon_gl = 0.0000018e-3 # S/cm2
axon_gna = 0.300e-3 # S/cm2
axon_gkdr = 0.0525e-3 # S/cm2
amp = 0 # amplitude in nA
dur = 600 # duration in ms
delay = 200 # delay in ms


activemodel(tstop,soma_gl,gcas,gcat,gkA,gkca,soma_gkdr,taucas,axon_gl,axon_gna,axon_gkdr,amp,dur,delay)


Our logic is based on there theoretically being two split nodes, one connecting the end of the soma to the axon1, and two connecting the end of axon1 to the beginning of axon2. However, we are unsure of how we should implement multisplit with respect to the splitnodes.

Would help if someone can help with the logic in the code above or point us to simple implementations of this type to parallelize single cell sections.

If you would like to reach out:
Atul Krishnadas: akyvv@unsystem.edu
Matthew Carroll: mjc6r9@umsystem.edu
ted
Site Admin
Posts: 6286
Joined: Wed May 18, 2005 4:50 pm
Location: Yale University School of Medicine
Contact:

Re: Parallelizing sections within a single cell in python without using GUI

Post by ted »

"As proof of concept" first create a non-parallelized implementation of your model cell, and work with that a bit to make sure that (1) parameter values are correct, and (2) simulation results are what you expected to observe.
soma.cm = 9e-3 # membrane capacitance uF/cm2
Experimental evidence for this value? AFAIK experimentally-based estimates of specific membrane capacitance range from 0.6 to about 2 uf/cm2.
soma.Ra = 10000 # ohm-cm
Cytoplasmic resistivity is harder to determine from experiments, but empirical observations and modeling studies that are closely linked to experimental results indicate values that range from ~ 30 ohm cm (35.4 in squid axon) to 200 - 250 ohm cm at the high end (most models of neurons from vertebrates assume values in the range 100 - 160 ohm cm).
h.nrnmpi_init()
pc = h.ParallelContext() . . .
Is parallelization even necessary, and if so, what style of parallelization is most appropriate--thread-parallelization, bulletin-board parallelization, or multisplit? Programmer time is more valuable than CPU time, and anything that complicates implementation or debugging will cost human time and effort.

Multithreaded simulation is easiest to implement--requires no change to the structure of the model setup code--and allows interactive simulation, which is extremely useful for debugging and exploration of model behavior (a key preliminary to e.g. optimization or parameter space exploration). The tradeoff is that speedup is sublinear, i.e. execution on N+1 cores may be slower than (N+1)/N times faster than execution on N cores, and even for very large models the maximum useful value of N is about 15.

Bulletin-board parallelization is best for embarassingly parallel problems (you want to run many simulations with the same model, e.g. to explore parameter space or quantify effects of stochasticity). It has much less overhead than multithreaded execution, so speedup tends to be linear with N as long as all processors are busy. It requires only a small amount of additional programming effort. However, you lose interactivity.

Multisplit simulation of individual cells is best suited to projects that involve long simulations (tstop is many thousands of ms or longer) of very complex model cells. This requires significant programming effort, and of course there is no interactivity.

A reasonable strategy would be to first see if you can get reasonable run times with single thread or multithreaded simulation. If not, then consider whether bulletin-board parallel simulation will do what you need, before moving on to multisplit parallelization of an individual cell.
bll5z6
Posts: 29
Joined: Thu May 26, 2016 10:27 am

Re: Parallelizing sections within a single cell in python without using GUI

Post by bll5z6 »

@ted, is there some documentation or an example of how to implement multithreaded simulation without the NEURON GUI?
ted
Site Admin
Posts: 6286
Joined: Wed May 18, 2005 4:50 pm
Location: Yale University School of Medicine
Contact:

Re: Parallelizing sections within a single cell in python without using GUI

Post by ted »

I think this reply might be useful to both of you.
bll5z6 wrote: Sun Apr 03, 2022 9:51 amis there some documentation or an example of how to implement multithreaded simulation without the NEURON GUI?
No. I am aware of one published article that implemented multithreaded simulation without mapping NEURON's "parallel computing" GUI interface (NEURON Main Menu / Tools / Parallel Computing) to the screen--see ModelDB entry 144054, and in particular mctavish_syncbylocation/src/parinit.hoc. What it does is beautiful (compile the mod files, then execute mosinit.hoc and follow the instructions). That said, the code was written 10 years ago, contains some additional "runtime reporting" statements that you may not need, and I am not sure that it follows current best practices.

Fortunately
1. NEURON's GUI tools are just interfaces to its computational engine and data structures, which exist regardless of whether you're using the GUI or not. The GUI is just a convenient way of setting switches and configuring those data structures without you having to write any code. When you save a configured GUI tool to a .ses file, that file contains not only what is necessary to redraw the tool on the screen, it also contains sufficient information to properly configure NEURON's computational engine and recreate the data structures that do the actual work.
2. NEURON's GUI tools will work even if your own model specification is written in Python. So if you want to use multithreaded simulation, instead of
from neuron import h
just do
from neuron import h, gui
and then you will be able to bring up a parallel computing GUI tool and a RunControl panel and use those interactively to (A) discover how many usable cores are on your PC or Mac, (B) see how much speedup you get from efficient cache usage alone, and (C) discover the optimal number of threads to use for best load balance and speedup
3. You can then save your optimally configured parallel computing GUI tool to a .ses file (look in the FAQ list at www.neuron.yale.edu -- link to FAQ will be on the Documentation page, and the particular FAQ item of interest is called I've used the NEURON Main Menu to construct and manage models. How can I save what I have done? -- to see how to save selected GUI panels to a .ses file). The code in that .ses file will contain not only the info that NEURON needs to recreate that GUI tool on your desktop, but also to recreate the internal data structures that you configured when you used the GUI tool in the first place. And you can make that happen just by using h.load_file() (or load_file() if you're working with hoc) to execute the contents of that .ses file.

Which brings My suggestion:
Work your way through this short exercise from the NEURON Course https://neuron.yale.edu/neuron/docs/mul ... lelization, which introduces multithreaded computing and shows you how to use the parallel computing tool to discover how many usable cores are available on your PC or Mac,

Then implement your model in Python but be sure to
from neuron import h, gui
so the NEURON Main Menu toolbar is present.

After your model specification etc. is complete, bring up a RunControl panel and a "Parallel Computing" GUI tool. Run some interactive simulations to discover whether cache efficient helps, and how many threads you want to use. Save the configured parallel computing tool to a .ses file; you might want to call this file something like threads.ses. Finally, in your own code AFTER model specification is complete but before the code that launches simulations, insert
load_file("threads.ses") // if you're using hoc
or
load_file("threads.ses") # if you're using python

And you should be good to go.
Post Reply