Here is an example used for computing correlations.

It's almost a direct translation of several hundred

lines of Fortran 77.

If you can bare to look, please help, even it is simply

to say "yo bonehead, why didn't you just require 'some_library'?"

Thanks,

## ···

--

Bil

http://fun3d.larc.nasa.gov

cat << EOF > correlation_test.rb

require 'test/unit'

require 'correlation'

class TestCorrelation < Test::Unit::TestCase

def test_1x1_perfectly_correlated

x = [ [0,1] ]

y = [ [0,1] ]

c, c2 = correlation(x,y)

assert_equal( [ [ 1.0 ] ], c )

assert_equal( [ [ 1.0 ] ], c2 )

end

def test_1x1_perfectly_anti_correlated

x = [ [0,1] ]

y = [ [1,0] ]

c, c2 = correlation(x,y)

assert_equal( [ [ -1.0 ] ], c )

assert_equal( [ [ 1.0 ] ], c2 )

end

def test_1x1_perfectly_uncorrelated

x = [ [0,1] ]

y = [ [0,0] ]

c, c2 = correlation(x,y)

assert_equal( [ [ 0.0 ] ], c )

assert_equal( [ [ 0.0 ] ], c2 )

end

def test_2x1_perfectly_correlated

x = [ [0,1], [0,1] ]

y = [ [0,1] ]

c, c2 = correlation(x,y)

assert_equal( [ [ 1.0 ], [ 1.0 ] ], c )

assert_equal( [ [ 0.5 ], [ 0.5 ] ], c2 )

end

def test_2x1_perfectly_correlated_and_uncorrelated

x = [ [0,1], [0,0] ]

y = [ [0,1] ]

c, c2 = correlation(x,y)

assert_equal( [ [ 1.0 ], [ 0.0 ] ], c )

assert_equal( [ [ 1.0 ], [ 0.0 ] ], c2 )

end

def test_1x2_perfectly_anti_correlated

x = [ [0,1] ]

y = [ [1,0], [1,0] ]

c, c2 = correlation(x,y)

assert_equal( [ [ -1.0, -1.0 ] ], c )

assert_equal( [ [ 1.0, 1.0 ] ], c2 )

end

def test_1x2_perfectly_correlated_and_uncorrelated

x = [ [0,1] ]

y = [ [0,1], [0,0] ]

c, c2 = correlation(x,y)

assert_equal( [ [ 1.0, 0.0 ] ], c )

assert_equal( [ [ 1.0, 0.0 ] ], c2 )

end

end

EOF

cat << EOF > correlation.rb

# Returns an array with correlation and y-normalized correlation^2

# Assumes two nested arrays of variables and samples. For example,

# here's a case with 3 samples of 3 input variables and 2 output variables:

#

# x = [ [a1,a2,a3], [b1,b2,b3], [c1,c2,c3] ]

# y = [ [r1,r2,r3], [s1,s2,s3] ]

def correlation(x,y)

fail 'size arrays unequal' if x.first.size != y.first.size # lame

n = x.first.size

sum_x = Array.new(x.size){0.0}

sum_x_sq = Array.new(x.size){0.0}

for i in 0..x.size-1

sum_x[i] = x[i].inject(0) { |sum, e| sum + e }

sum_x_sq[i] = x[i].inject(0) { |sum, e| sum + e**2 }

end

sum_y = Array.new(y.size){0.0}

sum_y_sq = Array.new(y.size){0.0}

for i in 0..y.size-1

sum_y[i] = y[i].inject(0) { |sum, e| sum + e }

sum_y_sq[i] = y[i].inject(0) { |sum, e| sum + e**2 }

end

sum_xy = Array.new(x.size){ Array.new(y.size){0.0} }

for k in 0..n-1

for i in 0..x.size-1

for j in 0..y.size-1

sum_xy[i][j] += x[i][k]*y[j][k]

end

end

end

corr = Array.new(x.size){ Array.new(y.size){0.0} }

for i in 0..x.size-1

for j in 0..y.size-1

dx = n*sum_x_sq[i] - sum_x[i]**2

dy = n*sum_y_sq[j] - sum_y[j]**2

corr[i][j] = ( n*sum_xy[i][j] - sum_x[i]*sum_y[j] ) / Math.sqrt(dx*dy) \

unless dx*dy==0.0

end

end

sum_corr_sq_y = Array.new(y.size){0.0}

for j in 0..y.size-1

for i in 0..x.size-1

sum_corr_sq_y[j] += corr[i][j]**2

end

end

corr_sq = Array.new(x.size){ Array.new(y.size){0.0} }

for i in 0..x.size-1

for j in 0..y.size-1

corr_sq[i][j] = corr[i][j]**2/sum_corr_sq_y[j] \

unless sum_corr_sq_y[j]==0.0

end

end

[corr, corr_sq]

end

EOF