In [5]:

    

import numpy as np
import pandas as ps

NUM_MATRICES = 1000

def rand_ortho_matrix(eigen_values):
    size = len(eigen_values)
    A = np.mat(np.random.random((size,size)))
    Q, R = np.linalg.qr(A)
    
    return Q * np.diag(eigen_values) * Q.T


    

In [2]:

    

def matrix_correlation(matrix1, matrix2):
    tr = matrix1.shape[0]
    correlation = np.corrcoef(matrix1[np.triu_indices(tr,1)], matrix2[np.triu_indices(tr,1)])[1, 0]
    return correlation

def random_skewers(matrix1, matrix2, num_vectors=1000):
    traits = matrix1.shape[0]
    rand_vec = np.random.multivariate_normal(np.zeros(traits),
                                             np.identity(traits, float),
                                             num_vectors).T
    delta_z1 = np.dot(matrix1, rand_vec)
    delta_z2 = np.dot(matrix2, rand_vec)

    ndelta_z1 = delta_z1/np.sqrt((delta_z1*delta_z1).sum(0))
    ndelta_z2 = delta_z2/np.sqrt((delta_z2*delta_z2).sum(0))

    return np.mean(np.diag(np.dot(ndelta_z1.T, ndelta_z2)))


    

In [3]:

    

rand_eigen = np.array(10*np.random.random(n_traits))

ms = []

for i in xrange(1000):
    ms.append(rand_ortho_matrix(rand_eigen))


    

In [6]:

    

rs0 = []

for i in xrange(NUM_MATRICES):
    rs0.append(random_skewers(ms[0], ms[i]))


    

    




---------------------------------------------------------------------------
ValueError                                Traceback (most recent call last)
<ipython-input-6-0ae2dbb1731c> in <module>()
      2 
      3 for i in xrange(NUM_MATRICES):
----> 4     rs0.append(random_skewers(ms[0], ms[i]))

<ipython-input-2-2902de40a553> in random_skewers(matrix1, matrix2, num_vectors)
     12     delta_z2 = np.dot(matrix2, rand_vec)
     13 
---> 14     ndelta_z1 = delta_z1/np.sqrt((delta_z1*delta_z1).sum(0))
     15     ndelta_z2 = delta_z2/np.sqrt((delta_z2*delta_z2).sum(0))
     16 

/home/walrus/.virtualenvs/py/local/lib/python2.7/site-packages/numpy/matrixlib/defmatrix.pyc in __mul__(self, other)
    328         if isinstance(other,(N.ndarray, list, tuple)) :
    329             # This promotes 1-D vectors to row vectors
--> 330             return N.dot(self, asmatrix(other))
    331         if isscalar(other) or not hasattr(other, '__rmul__') :
    332             return N.dot(self, other)

ValueError: objects are not aligned

In [23]:

    

rand_pop = np.random.multivariate_normal(np.zeros(n_traits), t_Ne*p, n_pop)


    

In [28]:

    

jklnp.cov(rand_pop, rowvar=0)


    

    
Out[28]:





(10, 10)

In [24]:

    

rand_pop.shape


    

    
Out[24]:





(100, 10)

In [ ]: