nikkip wrote:I'm using CVODE for the first time
First a question for you: do you need cvode.record? No need for that unless you're using the local variable time step method, and you won't be using that unless (1) your model has multiple cells, and (2) you have found that simulations run quicker or are significantly more accurate when executed with the local variable time step method.
The following answers assume that you do not need cvode.record.
I'm running my simulation multiple times (parameter sweep). Should I call cvode.record(...) just once (outside my simulation procedure) or for every simulation (within the simulation procedure)?
I don't know what you mean by "simulation procedure." In the context of computational modeling it is most useful to think of the following kinds of code:
model setup code--specifies the properties of the physical system that will be represented in the computer
instrumentation code--specifies signal sources or inputs to the physical system, and which signals or variables are to be observed or recorded
simulation control code--specifies how simulations will be initialized, and what happens in the course of simulation execution (duration of simulation, integration method, dt or error tolerance, conditions or times at which perturbations will occur (think "virtual experimental protocol"))
In general, properly written model setup and instrumentation code only needs to be executed once. This is certainly the case for statements that specify which variables are to be recorded. That said, if you are going to be writing code that contains unfamiliar features of NEURON, it is a very good idea to first try those features in the context of a toy problem (a simple task involving a very simple model for a very short run time) just to make sure that there will be a close match between what you want and what your code actually does.
2) I want to use cvode.record(...) to record the current of my IClamp point process.
Unless you're using local variable timestep integration, use the Vector class's record method. In either case, make a toy model (e.g. a single compartment with an IClamp attached to it), and see if the values recorded are what you expected based on the IClamp's parameters. The documentation of Vector.record provides several examples, the first of which is
vdest.record(&var)
so for anything you might want to record, first ask "what is the name of the variable in hoc, and what is the name of the Vector to which I want to record it?" and that should give you a good idea of what the necessary statement is. Examples:
Suppose you have an IClamp called stim and you want to record its time course to a vector called ivec. The current that stim delivers is called stim.____ so the statement is ivec.record(&_______).
Suppose you also want to record time at each step in the simulation to a vector called tvec (a good idea, because you can then use the Vector class's plot method to show the time course of stim.i in a graph. Time is called t, so the statement is tvec.record(&______).
Now make a toy model and try out these statements.
3) I would also like to use cvode.record(...) to record Vm.
The documentation of Vector.record tells what to do.
Actually, that's not a vector--it's an example of the use of array notation to manage multiple things (sections in this case) that have the same base name. The complete name of any one of your nodes is node
where i is either a whole number or a variable whose value lies in the range 0..axonnodes-1. So if the way to record v at the midpoint of a section called terminal to a vector called dv is
dv.record(&terminal.v(0.5))
what do you think would be the way to record terminal.v(0.5) if the section's name is node[3] instead of terminal? How would you test this idea?
With regard to the examples you provide, the problem is a very understandable misinterpretation of the documentation. Referring to the documentation for Vector.record, here are the key points:
vdest.record(&var) means that the values of var will be captured to vdest after finitialize() and each fadvance().
vdest.record(&var, Dt) means that values will be captured to vdest at the first simulation times that are >= i*Dt where i = 0, 1, 2 . . .. If Dt is a multiple of dt, this will capture every "Dt/dt"th value generated in the simulation, allowing for roundoff error here or there.
vdest.record(&var, tvec), where tvec is monotonically increasing, means that values are captured to vdest at the first simulation times that are >= the values of the elements in tvec. In other words, tvec is a user-supplied vector of times at which you want the value of var to be sampled, again allowing for roundoff error here and there.
So your code examples don't capture anything because the vector of sample times is empty. But why bother with that syntax, unless you are using adaptive integration and really want to capture values at specific times? And even then you only need to use Vector record unless the model has multiple cells and the simulation is being executed with local variable time steps.
4) I’m using extracellular current point source stimulation (phi = I/(4*pi*sigma*r) at each node, where I is a current pulse from an IClamp). I want to set up the extracellular potential time course for each node and use myvec.play(…) to apply them. I’ve tried looking at the folder "extracellular_stim_and_rec", but I’m having trouble grasping the pieces that I need. From calcrxc.hoc, I understand that I can create a vector (rx_xtra) that contains the transfer impedance from the extracellular point source to each node of my axon. Then, using rx_xtra, I could just loop through all the nodes and all time points with a fixed time step…
Code: Select all
for K = 0, axonnodes-1 {
for J = 0, ntsteps – 1 {
// compute Ve at each node at each time point
if ((J*dt > delay) && (J*dt < delay + dur)) {
Ve[K].x[J] = I/(4*pi*sigma*r.x[K])
}
}
}
Is there a more efficient way, in which I’d actually use the IClamp rather than figuring out if the stimulation is on or not manually as I did above?
Then I think I’d have to apply (i.e. play) the potentials with something like…
Code: Select all
for K = 0, axonnodes-1 {
Ve[K].play(&node[K].e_extracellular(0.5), dt)
}
I will also use cvode.event(t) to add events signaling the step up and step down of my current pulse.
This is about a very different topic, so I will split it off into a different thread and answer it there at a later time.