Page 1 of 1
Data analysis
Posted: Fri Feb 24, 2012 9:17 pm
by erhervert
Hello again,
I'd like to know if there is a way to customize the output ASCII file obtained from a simulation for further analysis. For example, if I execute my model n times and then I used the print & file window manager to select the V(t) plot, what I get is a file with two columns, one x and one y containing all the results. But, what I'd like to do is to get a file with a single x (time column) and several y's (y0, y1, y2,..., yn), one y from each run. Do you know if it's possible to automate that? Or if there is a way to average all the vectors?
Thanks in advance
Re: Data analysis
Posted: Sat Feb 25, 2012 12:56 am
by ted
Set up Vector recording of the variable you are interested in.
Create another Vector that will eventually hold the average over several runs.
For the sake of discussion, let's say the two Vectors are called tmp and avg.
Now automate the process of generating several runs and calculating the average.
Code: Select all
proc batchavg() { local i
avg = new Vector() // discard whatever is already in avg
for i = 1,$1 {
run()
if (i==1) {
avg = tmp.c // copy tmp to avg
} else {
avg.add(tmp) // add tmp to avg, element by element
}
avg.div($1) // divide each element by number of runs to get the average
}
}
Example of usage:
batchavg(5)
runs 5 simulations and fills a Vector called avg with the average time course of the recorded variable.
Re: Data analysis
Posted: Sat Feb 25, 2012 11:51 am
by erhervert
Ok, I incorporate that into my program and this is what I got:
Code: Select all
///////////////////////////
//* simulation control *//
///////////////////////////
dt = 0.025
tstop = 50
v_init = -65
objref r
r = new Random()
r.normal(.5, .05)
objref avg
proc batchavg() { local i, index
for i=0, $1 {
index = r.repick()
syn.onset = (((10-5)*index)+5)
run()
}
avg = new Vector()
if (i==1) {
tmp = avg.c()
} else {
avg.add(tmp)
}
avg.div($1)
}
Apparently it is working well, but I still having the same two colums (x,y). Do I have to create a second graph to plot the vector avg vs time to see the resulting data? How can I do it? I have been trying to plot the vector as follows, based on the programmers reference example (
http://www.neuron.yale.edu/neuron/stati ... .html#plot) but it is not working:
Code: Select all
objref avg, g
g = new Graph()
g.size(0,50,-80,200)
avg.plot(g, 0.025)
graphList[0].append(g)
I'm not sure if it is correct.
Thanks
Re: Data analysis
Posted: Sun Feb 26, 2012 10:26 am
by ted
First a comment: if you don't prevent negative values from being assigned to your synapse's start time, you're going to have some runs that produce nonuseful results.
If you want to see the averaged result in a graph, you have to plot it _after_ all of the runs have finished--that is, after batchavg is executed. The graphLists are irrelevant to this. Use the Vector class's plot method. Also, you will probably want the x values to be the times that correspond to the data points, so it is necessary to record t to its own Vector (do that with code next to the statements that set up recording of v or whatever other variable you are recording).
The same considerations apply if you want to see a printout of time and the averaged values.
Re: Data analysis
Posted: Sun Feb 26, 2012 5:56 pm
by erhervert
First a comment: if you don't prevent negative values from being assigned to your synapse's start time, you're going to have some runs that produce nonuseful results.
Sorry, I did it already using the tsynon() function as you recommend but I posted the uncorrected code (by mistake), this is the good one:
Code: Select all
objref r
r = new Random()
{r.normal(.5, .05)}
func tsynon() {
ton = -1
while (ton<0) {
ton = (((10-5)*r.repick())+5)
}
return ton
}
I also did the plot of average vs time using the vector class's plot method. In addition, I incorporate some code to write the results in a txt file at the end of the simulation. But I'm not sure if it is calculating the average correctly, the resulting signal is smaller in size than a single AP and I was expecting something similar to a compound action potential (the sum of the action potentials).
This is the code, could you tell me what do you think?:
Code: Select all
objref avg, time, savdata, tempmatrix, g
time = new Vector()
{time.record(&t)}
avg = new Vector()
{avg.record(&soma.v(0.5))}
proc batchavg() { local i
for i=0, $1 {
syn.onset = tsynon()
run()
}
if (i==1) {
tmp = avg.c()
} else {
avg.add(tmp)
}
g = new Graph()
g.size(0,50,-10,5)
avg.plot(g, time)
avg.div($1)
savdata = new File()
savdata.wopen("results.dat")
savdata.printf("t avg\n")
tempmatrix = new Matrix()
tempmatrix.resize(avg.size(),2)
tempmatrix.setcol(0, time)
tempmatrix.setcol(1, avg)
tempmatrix.fprint(savdata, " %g")
savdata.close()
}
I'd like to know also if there is a way to make measures on the average trace. I'm interested in things like the peak amplitude, dv/dt and the area under the curve.
Thanks for your help Ted.
Re: Data analysis
Posted: Mon Feb 27, 2012 9:36 am
by ted
Looks good. My only suggestion would be to split the "write data to file" stuff out of batchavg and into its own proc. Then you'd pull them together with
Code: Select all
strdef fname
fname = "name_for_your_file.dat"
proc doall() {
batchavg(...)
savetofile(fname)
}
Rationale: modular code is easier to write and debug. Also makes it easier to reuse for other purposes.
The Vector class has lots of procs and funcs that are very useful for analyzing data. Check it out in the Programmer's Reference (the "Quick Index" is particularly convenient).
Re: Data analysis
Posted: Mon Feb 27, 2012 2:29 pm
by erhervert
Done
Code: Select all
objref savdata, tempmatrix
proc savetofile() {
savdata = new File()
savdata.wopen("results.dat")
savdata.printf("t avg\n")
tempmatrix = new Matrix()
tempmatrix.resize(avg.size(),2)
tempmatrix.setcol(0, time)
tempmatrix.setcol(1, avg)
tempmatrix.fprint(savdata, " %g")
savdata.close()
}
proc doall() {
batchavg(10)
savetofile()
}
Thanks again Ted.
Re: Data analysis
Posted: Mon Feb 27, 2012 9:51 pm
by ted
Hey, I've got it easy. You're the one who has to actually write code and fix it when it doesn't work.
Re: Data analysis
Posted: Mon Feb 27, 2012 11:29 pm
by erhervert
I found it addictive, feels great when the things start working good.