Re: [Matplotlib-users] matplotlib.mlab.corrcoef() not working?

[John H. <jdh...@ac...>]

>>>>> "Sebastian" == Sebastian Stark <st...@tu...> writes:

    Sebastian> How is corrcoef() supposed to work? When I do the
    Sebastian> following:

    Sebastian> --------------------------------------- from numarray
    Sebastian> import * from numarray import random_array from
    Sebastian> matplotlib.mlab import corrcoef randn =
    Sebastian> random_array.standard_normal

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

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