OpenEuphoria

Re: Correlation Coefficient

• Posted by Juergen Luethje <j.lue at gmx.d?> Jul 29, 2007
• 538 views

I tested the speed of 4 different functions for calculating the
standard deviation (interpreted by exw.exe, Eu 3.0.2):

without type_check

function sum (sequence data)
   atom t
     
   t = 0
   for i = 1 to length(data) do
      t += data[i]
   end for 
   return t
end function

------------------------------------------------------------------------

function std_dev_1 (sequence data)
   -- from Rob's post
   atom avgx
    
   avgx = sum(data) / length(data)
   return sqrt( sum(power(data - avgx, 2)) / length(data) )
end function

function std_dev_2 (sequence data)
   atom avgx
   
   avgx = sum(data) / length(data)
   data -= avgx
   data *= data
   return sqrt( sum(data) / length(data) )
end function

function std_dev_3 (sequence data)
   -- Note: _Two_ loops are required here, one loop in
   --       sum(data) and another one in this function.
   atom avgx, d, sum_d
   
   avgx = sum(data) / length(data)
   sum_d = 0
   for i = 1 to length(data) do
      d = data[i] - avgx
      sum_d += d*d
   end for
   return sqrt( sum_d / length(data) )
end function

function std_dev_4 (sequence data)
   -- Note: This calculation does not use the average,
   --       so it's not necessary to calculate sum(data)
   --       first, hence only _one_ loop is required here!
   atom x, sum_q, sum_x
   
   sum_q = 0
   sum_x = 0
   for i = 1 to length(data) do
      x = data[i]
      sum_q += x*x
      sum_x += x
   end for
   return sqrt( (sum_q - sum_x*sum_x/length(data)) / length(data) )
end function

------------------------------------------------------------------------

sequence a
atom t1, t2, t3, t4

a = rand(repeat(100, 10000000))

t1 = time()
? std_dev_1(a)
t1 = time() - t1

t2 = time()
? std_dev_2(a)
t2 = time() - t2

t3 = time()
? std_dev_3(a)
t3 = time() - t3

t4 = time()
? std_dev_4(a)
t4 = time() - t4

puts(1, "\n")
               -- Results on my PC: 
? t1           -- 5.87 
? t2           -- 1.55 
? t3           -- 1.06 
? t4           -- 0.55 

Regards,
   Juergen