How is corrcoef() supposed to work? When I do the

> following:

> --------------------------------------- from numarray

> import * from numarray import random_array from

> matplotlib.mlab import corrcoef randn =

> random_array.standard_normal

> x = randn((10,4)) r = corrcoef(x)

> ---------------------------------------

Hmm, I don't know how this slipped in untested. corrcoef definitely

works for corrcoeff(x,y), but is obviously broken for matrices.

Here is a fixed implementation, tested against matlab for MxN matrices

with 5 significant digits of agreement in the correlation matrix

output with Numeric and numarray.

Thanks for catching it.

As you noted, there was a non-printable character n the buffer, which

crept in when I borrowed some code from Fernando Perez. Damned

Colombians and their accented names!

Let me know if you encounter further problems!

JDH

def corrcoef(*args):

"""

corrcoef(X) where X is a matrix returns a matrix of correlation

coefficients for each numrows observations and numcols variables.

corrcoef(x,y) where x and y are vectors returns the matrix or

correlation coefficients for x and y.

Numeric arrays can be real or complex

The correlation matrix is defined from the covariance matrix C as

r(i,j) = C[i,j] / sqrt(C[i,i]*C[j,j])

"""

if len(args)==2:

X = transpose(array([args[0]]+[args[1]]))

elif len(args)==1:

X = args[0]

else:

raise RuntimeError, 'Only expecting 1 or 2 arguments'

C = cov(X)

if len(args)==2:

d = resize(diagonal(C), (2,1))

denom = sqrt(matrixmultiply(d,transpose(d)))

else:

dc = diagonal(C)

N = len(dc)

shape = N,N

denom = sqrt(matrixmultiply(diagonal_matrix(dc),

resize(dc, shape)))

r = divide(C,denom)

try: return r.real

except AttributeError: return r