plotting data as the simulation runs

[porio]

Hi,

Whenever I create a Graph from a .hoc code, I must explicitely add the data to be plotted and this is usually done after the simulation has run, using a vector that recorded the data of interest.
I would like to create a Graph that is updated during the simulation, like the graphs that are created using the GUI.

I dived into the GUI's hoc files and found the newPlotV() procedure. I can call it from my hoc and it does what I want, but then I tried to reproduce the procedure (using different object names but the same commands and looking recursively into the procedures), and it didn't work; the graphs (or plots?) that I create always have to be updated explicitely after the simulation.

Another question: if I call the newPlotV() command from my code (to use an easy solution), what would be the most convenient way of referring to the plot that is created in order to change the scale, add more data or change the plot color?

Thanks!

[Raj]

You can simply use the
newPlotV, newPlotS or newPlotI functions for plotting Voltage, State Variables and Currents respectively.

Make sure you have the standard GUI hoc loaded.

All three functions mentioned above store a reference to the graph in graphItem, so for example state variables can be plotted using newPlotS.

Code: Select all

newPlotS() 
objref MyGraph
MyGraph=graphItem

This code snippet generated a state plot and saves the reference to it in MyGraph.

Now you can add new variables to MyGraph using addexpr:

Code: Select all

MyGraph.addexpr("synapse.x", 2, 1)
MyGraph.addexpr("synapse.s", 3, 1)
		

Addexpr gives you the opportunity to plot extra variables and you can specify color and line style with it. Addexpr is a method of the Graph class and indeed graphItem is a reference to a Graph class object. The programmers reference on the neuron website gives a full description of this class so you should be able to find the information you are looking for there.

The difference between newPlotS, newPlotV, newPlotI is due to the numerical algorythm which uses a staggered time step, using these functions makes your program plot them at the time when they are known at the highest accuracy.

[ted]

The difference between newPlotS, newPlotV, newPlotI is due to the numerical algorythm which uses a staggered time step, using these functions makes your program plot them at the time when they are known at the highest accuracy.

The story is a bit more complicated.

Staggered time steps are used only by NEURON's Crank-Nicholson integrator, as part
of a clever scheme for achieving second order accuracy in time while simulating models
that involve HH-sytle mechanisms, without having to resort to iteration. This scheme
computes currents, voltages, and states, all to second order accuracy, BUT at slightly
different times--t-dt/2, t, and t+dt/2, respectively. You may have noticed that the GUI's
Voltage axis, Current axis, and State axis graphs plot variables with these temporal
offsets, when you're using fixed time step integration. Now you know why.

The kid in the front row with the propeller beanie and "Nerds rule" T-shirt pops up
his hand to object. "Doesn't this mess up plots generated with NEURON's default
integrator?"

Not really. The default integrator is an implicit Euler method, which is first order
accurate in time. Given a first order accurate solution computed over the interval
[t, t+dt], there is no "time of maximum accuracy" in this interval. It's just as accurate
(or inaccurate) to plot things shifted by half a time step, as to plot them all at the same
time.

When adaptive integration is used ("Use variable dt"), the staggered times go away--
currents, voltages, and states are all computed to the same degree of accuracy at
the same instant, and all three types of graphs automatically switch to plotting points
at the same instant.

All of this raises a big question for people who write hoc code to bring up graphs:
how make sure that their code plots things at the proper times? There's a very quick
and easy trick for doing this, which I'll describe in a separate posting to this thread.

[ted]

In a separate posting, I discussed NEURON's use of staggered time steps and touched
on what this implies for plotting variables.

All of this raises a big question for people who write hoc code to bring up graphs:
how make sure that their code plots things at the proper times? There's a very quick
and easy trick for doing this, which I'll describe in a separate posting to this thread.

Here's the trick: use the GUI.

Right. Use NEURON Main Menu / Graph / Voltage axis
(or Current axis, or State axis, as appropriate) to bring up a graph of the desired type.
Then use that graph's Plot what? (and, if necessary, Delete) to add whatever variable(s)
you need to the graph--taking care to plot voltages on a Voltage axis graph, currents
on a Current axis graph, and States on a State axis graph. Save the configured graph
to a session file. When you need it, read it in with a load_file statement. You can also
edit a session file, and mine reusable code from it.

"What about ionic conductances?" asks propellerhead.

Good question. If you have a mechanism with a voltage- or ligand-gated conductance,
the conductance will be an ASSIGNED variable, just like currents are, so you just
plot it in a Current axis graph.

For example, suppose we have a model with a sodium current in its default section.
We bring up a Current axis graph, use its Plot what? to specify ina, and save the
graph to a session file. Here's the key part of the session file:

Code: Select all

{
save_window_ = new Graph(0)
save_window_.size(0,5,-1,1)
scene_vector_[5] = save_window_
{save_window_.view(0, -1, 5, 2, 523, 375, 300.48, 200.32)}
graphList[1].append(save_window_)
save_window_.save_name("graphList[1].")
save_window_.addexpr("ina", 1, 1, 0.8, 0.9, 2)
}

It would be trivial to turn this into a graph of ik--just change the first argument to the
addexpr statement. Or we could add ik to the plot list by inserting this statement:
save_window_.addexpr("ik", 1, 1, 0.8, 0.85, 3)
(changes in bold). Then we save it back to a session file and load it whenever
we need it.

Or we could just steal reusable code that, after some modification, we can insert
directly into our own hoc code. For example, this creates a new graph of ina:

Code: Select all

objref g
g = new Graph(0)
g.size(0,5,-1,1)
g.view(0, -1, 5, 2, 523, 375, 300.48, 200.32)}
graphList[1].append(g)
g.addexpr("ina", 1, 1, 0.8, 0.9, 2)

"That's an awful lot of code. And I don't like where this graph pops up, either."

First, we didn't have to write most of this code--we just llifted it from the ses file
and made a few minor changes. Second, we could make most of these changes
by analogy, i.e. without having to read the Programmer's Reference! Thrid, this code
gives us full control over graph location, size, and axis scaling. If you don't like where
it is or how it looks, drag it around, resize it, use its GUI to change colors and brushes,
and then save it to a session file. Steal the desired coordinates (and other changes)
from that file, and use them to revise your hoc code.

[Raj]

I spent a few years of my life writing GUI code for a MS-windows application. I think the amount of code is actually quite modest, because the Graph class is nicely tailored to do the job.

So except for some subtleties related to different integration scheme's, about which (when needed) I'll try to write more carefully in the future, its a breeze.

When it comes to CVode I'm not sure whether this code will do the job in combination with a variable timestep. In my experience some variables just jump around like mad and even error messages result from having a graph open. This is probably due to settings like use_local_dt and/or condition_order(2) which make time local and no longer make time strictly incremental.

Using CVode's record method and plotting the resulting vectors after the simulations seems a more reliable approach in these situations.

Therefore some example code. You'll have to add some current(s) to this example to see something happening.

Code: Select all

create soma
access soma

objref cvode
cvode = new CVode()
	
IsActive=cvode.active(1)
fprint( "IsActive: %d \n", IsActive)

UseLocalDtOn=cvode.use_local_dt(1)
fprint( "Use Local dt: %d \n", UseLocalDtOn)
	
cvode.condition_order(2)
ConditionOrder=cvode.condition_order()
fprint( "Condition Order: %d \n", ConditionOrder)

objref yvec, tvec
yvec= new Vector()
tvec = new Vector()
cvode.record(&v(0.5), yvec, tvec)

objref g
g = new Graph(0)

xpanel("Plot Traces")
	 xbutton("Plot Soma Membrane Potential","yvec.plot(g,tvec)")
xpanel()

Which is a bit more code, but I wanted it to be clear from the code that CVode is involved. If you call yvec.plot directly from your code it will most likely be executed before you do a run, in which case the vector is empty and hence your graph will stay empty aswell. Therefore I added a panel from which you can plot the variable after executing (part of ) your simulation. I assume there is a way of calling these functions automatically at the end of a run, which would make life a lot easier. I did not dive into that issue yet, maybe later.

[ted]

When it comes to CVode I'm not sure whether this code will do the job in combination with a variable timestep. In my experience some variables just jump around like mad and even error messages result from having a graph open.

That's very strange; I haven't ever seen such behavior. Graphs spawned from the
GUI are supposed to work with all of NEURON's integration methods, without any
special intervention by the user. Saving such a graph to a session file actually writes
hoc code that will, when executed, recreate the graph--which should then behave
properly. Variables "jumping around" suggests that you may have uncovered a bug.
If you can provide some code that reproduces this effect, it would be helpful in tracking
the bug down and eliminating it.

This is probably due to settings like use_local_dt and/or condition_order(2) which make time local and no longer make time strictly incremental.

The standard run system is supposed to take care of graph updates properly, even
when local time steps are used and variables from different cells are plotted in the
same graph. I suppose it is possible to write graph plotting routines that bypass the
standard run system, and this might have the effects you describe.

An aside--
Several very productive modelers, who began using NEURON years ago, have
developed their own programming styles and code libraries for simulation control and
graph updating. You'll find many examples of this scattered around the web and in
ModelDB. These idiosyncratic, if effective, approaches work quite well in their hands,
and they have no need or motivation to do otherwise in the future. However, the vast
majority of NEURON users should take advantage of the power and flexibility of the
standard run system, which has the advantages of working seamlessly with all of
NEURON's integrators and GUI tools, and being very easy to customize. The standard
run system, and practical tips for using it to control simulations and plot results, are
discussed in chapters 6-8 of The NEURON Book.

An example of customizing the standard run system

I assume there is a way of calling these functions automatically at the
end of a run

There is, and it involves a minor customization of the standard run system.

1. Insert
load_file("nrngui.hoc")
at the top of your hoc code. This makes the standard run library available.
2. Insert

Code: Select all

proc run() {
  stdinit()
  continuerun(tstop)
  yvec.plot(g.tvec)
}

at the end of your hoc code (or at least after the load_file(nrngui.hoc) statement).
This replaces the standard run system's proc run() (which consists of just stdinit()
and continuerun(tstop)) with your own custom run().

Now if you type
run()
at the oc> prompt, the simulation will execute, and after it finishes, the captured data
will be plotted.

Two comments:
1. This is approach can be generalized for the purpose of automatically launching
postprocessing of any kind after the end of a simulation.
2. Adding a Run button to Raj's xpanel "is left as an exercise to the reader."