Statistic analysis using gnuplot (1)

Last time we dealt with the maximum, minimum and mean value of a statistic analysis problem. Today the standard deviation will be added.
The standard deviation is defined as the square root of the variance. With variance stand for the mean squared value minus the squared mean value. The mean value have been dealt last time. We just need to store it. Strange it is, although there is only one data point after the smooth, when I plot it to a data file using the table mode, I get two output points with the second one marked with "u" which means undefined value.(who can tell me why?) In a similar way we plot the mean of squared value to another data file. Next we read this two value out from the data file. And use them to calculate the standard deviation. Since there is only two data points in the file, we use f(x)=ax+b to fit the data, and the fitted result will be exact, then we can use f(x) to get these exact values. At last we plot the standard deviation on the graph. The following is our final plotting script.

reset
plot "rand_t.dat" u 1:2	#To get the max and min value
ymax=GPVAL_DATA_Y_MAX
ymin=GPVAL_DATA_Y_MIN
ylen=ymax-ymin
xmax=GPVAL_DATA_X_MAX
xmin=GPVAL_DATA_X_MIN
xlen=xmax-xmin
set table "mean.txt"	#put the mean value out
plot "rand_t.dat" u ((xmax+xmin)/2.0):($2) smooth unique w p
unset table
set table "mean_squared.txt"	#put the mean of squared value out
plot "rand_t.dat" u ((xmax+xmin)/2.0):($2**2) smooth unique w p
unset table
f(x)=a*x+b	#The fit functions
g(x)=c*x+d
fit f(x) "mean.txt" u 0:2 via a,b	#Read the mean and mean of squared value
fit g(x) "mean_squared.txt" u 0:2 via c,d
mean=f(0)	#mean value
mean_squared=g(0)	#mean of the squared value
standard_deviation=sqrt(mean_squared-mean**2)	#standard deviation
print "The mean value is ",mean		#print the mean and standard deviation
print "The standard deviation is ",standard_deviation
#plot
set term post eps enhanced color lw 1.5 font ",20"
set output "statistic.eps"
set xrange [xmin:xmax]
set yrange [ymin-0.5*ylen:ymax+0.5*ylen]
set xlabel "time"
set ylabel "Random Signal"
#The labels
set label 1 at (xmin+xmax)/2.,ymax "Maximum" offset 0,0.5
set label 2 at (xmin+xmax)/2.,ymin "Minimum" offset 0,-0.5
set label 3 at (xmin+xmax)/2.,mean "Mean" offset 0,0.5
set label 4 at (xmin+xmax)/2.,mean+3*standard_deviation \
             "Mean+3{/Symbol \163}" offset 0,0.5
set label 5 at (xmin+xmax)/2.,mean-3*standard_deviation \
             "Mean-3{/Symbol \163}" offset 0,-0.5
plot "rand_t.dat" u 1:2 w p pt 7 ps 0.5 notitle,\
     mean w l lt 2 notitle "Mean",\
     ymax w l lt 3 notitle "Maximum",\
     ymin w l lt 3 notitle,\
     mean+3*standard_deviation w l lt 4 notitle,\
     mean-3*standard_deviation w l lt 4 notitle
#the six plot stand for raw data, mean, maximum, minimum,
#3-sigma upper line and 3-sigma lower line

At last we get figure like follows.

Statistic analysis using gnuplot

Statistic analysis using gnuplot (0)

I will talk about statistic analysis using gnuplot in this article. The following contents are covered---maximum, minimum, mean value, standard deviation.

We begin with the maximum and minimum as they are the simplest. There are two gnuplot defined variables (to see all the gnuplot variables use command "show variable all"), GPVAL_DATA_Y_MAX and GPVAL_DATA_Y_MIN. They are the maximum and minmum y value in the data file which you just used to plot. So after plotting you look over the gnuplot defined variables and you find the maximum and minimum values.

The mean value is a bit difficult. There is no such a gnuplot defined value. To find it out we will play a trick. There is smooth option called "unique" which makes the data monotonic in x and points with the same x-value are replaced by a single point having the average y-value. So it is appropriate for finding the mean value. We use command

plot "data_t.dat" u (constant-value):($n) smooth unique w point

to plot the mean value. Here constant-value is an arbitrary constant and $n is the column which contains the random data. Note that there is only one point plotted using the above command, so plot style line can not be used. To plot the the mean value using a line we may need to use xerrobars plot style. (Next time I will give a method to plot it using plot style line.)

plot "data_t.dat" u (constant-value):($n):(xerrobar-length) smooth unique w xerrorbars

Now we come to plot the data points, maximum, minimum and mean value to a picture. And standard deviation is left next time.

reset
plot "rand_t.dat" u 1:2	#To get the max and min value
ymax=GPVAL_DATA_Y_MAX
ymin=GPVAL_DATA_Y_MIN
ylen=ymax-ymin
xmax=GPVAL_DATA_X_MAX
xmin=GPVAL_DATA_X_MIN
xlen=xmax-xmin
#plot
set term png
set output "statistic.png"
set xrange [xmin:xmax]
set yrange [ymin-0.5*ylen:ymax+0.5*ylen]
set xlabel "time(ms)"
set ylabel "Random Signal(Arbitary Unit)"
plot "rand_t.dat" u 1:2 w p pt 7 ps 0.5 notitle,\
     "rand_t.dat" u (xmax+0.1*xlen):($2):(1.1*xlen)\
     smooth unique w xerrorbars notitle,\
     ymax w l lt 3 notitle,\
     ymin w l lt 3 notitle
#plot raw data, mean value, maximun and minimun respectively
#There is only one data point for the mean value. To make it 
#looks like a line, we plot it using a xerrorbars plot style.

rand_t.dat can be downloaded here. Picture statistic.png is shown below.

Statistic analysis using gnuplot--maximum,minimum and mean value

Statistic analysis and histogram plotting using gnuplot

Given a data file containing a set of data, count how many datas locate in intervals [a1:a2],[a2:a3]... respectively, then plot the result into a histogram. This a common problem in statistics and exactly what we will do in this article.

Firstly, let us see how to map the data into intervals. There is a function "floor(x)" which return the largest integer not greater than its argument. So function floor(x/dx)*dx will map x into one of the intervals [-n*dx:-(n-1)*dx],[-(n-1)*dx:-(n-2)*dx]...[(n-1)*dx:n*dx].

Now we come to count the data number in each interval. In gnuplot there is a smooth option called "frequency". It makes the data monotonic in x. Points with the same x-value are replaced by a single point having the summed y-values. Using this property, we can count the data numbers in the intervals.

At last we plot our result using boxes plot style.

The main idea have introduced. It is time to write the plotting script.

reset
n=100 #number of intervals
max=3. #max value
min=-3. #min value
width=(max-min)/n #interval width
#function used to map a value to the intervals
hist(x,width)=width*floor(x/width)+width/2.0
set term png #output terminal and file
set output "histogram.png"
set xrange [min:max]
set yrange [0:]
#to put an empty boundary around the
#data inside an autoscaled graph.
set offset graph 0.05,0.05,0.05,0.0
set xtics min,(max-min)/5,max
set boxwidth width*0.9
set style fill solid 0.5 #fillstyle
set tics out nomirror
set xlabel "x"
set ylabel "Frequency"
#count and plot
plot "data.dat" u (hist($1,width)):(1.0) smooth freq w boxes lc rgb"green" notitle

We use a data file (download from here) which contains 10000 normally distributed random numbers and get a graph like the follow one.

statistic histogram plotting using gnuplot

Manipulate data using ternary operator in gunplot

In gnuplot there is an operator "? :" which is ternary. It behaves as it does in C language. For example, "a?b:c" means if a is true then return 'b', otherwise return 'c'. Using this operator we can do some complex operation on input data. Let us see an example.

We have a data file as follows:

#x  y1   y2
0.1	0.2	0.2
0.2	0.4	0.4
0.3	0.6	0.6
0.4	0.8	0.8
0.5	1	1
0.6	1.2	1.2
0.7	1.4	1.4
0.8	2.4	0.6
0.9	2.7	0.675
1	3	0.75
1.1	3.3	0.825
1.2	3.6	0.9
1.3	3.9	0.975
1.4	4.2	1.05
1.5	4.5	1.125

We want to plot the points y1=y2 with blue pints, and other points y1!=y2 (y1 not equals y2) with some other colors. To realize this goal, we can use ternary operator. According to the definition of ternary operator, expression "$2==$3?$2:1/0" will return the points satisfy $2=$3, and when $2!=$3 it returns undefied value(1/0 is an undefied value). So when we use command "plot 'data.dat' u 1:($2==$3?$2:1/0) w p", we plot the points with y1=y2. In a simillar manner, we can plot the points with y1!=y2. Now we come to our final plotting script.

reset
set term png font "Times,20"
set output "ternaty.png"
set xlabel "x"
set ylabel "y"
plot "data.dat" u 1:($2==$3?$2:1/0) w p lc rgb"blue" notitle,\
     "data.dat" u 1:($2!=$3?$2:1/0) w p lc rgb"green" notitle,\
     "data.dat" u 1:($2!=$3?$3:1/0) w p lc rgb"red" notitle

The picture comes out to be like follows.

Manipulate data using ternary operator in gnuplot