Data analysis

Anything that doesn't fit elsewhere.
erhervert
Posts: 13
Joined: Mon Nov 21, 2011 10:42 pm
Location: México City

Data analysis

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?

ted
Posts: 5704
Joined: Wed May 18, 2005 4:50 pm
Location: Yale University School of Medicine
Contact:

Re: Data analysis

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
for  i = 1,\$1 {
run()
if (i==1) {
avg = tmp.c // copy tmp to avg
} else {
}
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.

erhervert
Posts: 13
Joined: Mon Nov 21, 2011 10:42 pm
Location: México City

Re: Data analysis

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.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

ted
Posts: 5704
Joined: Wed May 18, 2005 4:50 pm
Location: Yale University School of Medicine
Contact:

Re: Data analysis

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.

erhervert
Posts: 13
Joined: Mon Nov 21, 2011 10:42 pm
Location: México City

Re: Data analysis

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)
}
}
``````
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 {
}
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.

ted
Posts: 5704
Joined: Wed May 18, 2005 4:50 pm
Location: Yale University School of Medicine
Contact:

Re: Data analysis

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).

erhervert
Posts: 13
Joined: Mon Nov 21, 2011 10:42 pm
Location: México City

Re: Data analysis

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.

ted
Posts: 5704
Joined: Wed May 18, 2005 4:50 pm
Location: Yale University School of Medicine
Contact:

Re: Data analysis

Hey, I've got it easy. You're the one who has to actually write code and fix it when it doesn't work.

erhervert
Posts: 13
Joined: Mon Nov 21, 2011 10:42 pm
Location: México City

Re: Data analysis

I found it addictive, feels great when the things start working good.