On Wed, Feb 15, 2012 at 10:52 PM, <josef.p...@gmail.com> wrote: > Doing a bit of browsing in the numpy tracker, I found this. From my > search this was not discussed on the mailing list. > > http://projects.scipy.org/numpy/ticket/1842 > > The multivariate normal random sample is not always the same, even > though a seed is specified. > > It seems to alternate randomly between two patterns, the sum of the > random numbers is always the same. > > Is there a randomness in the svd? Given the constant sum, I guess the > underlying random numbers are the same. > > (windows 7, python 32bit on 64 machine, numpy 1.5.1) > > import numpy > > d = 10 > alpha = 1 / d**0.5 > mu = numpy.ones(d) > R = alpha * numpy.ones((d, d)) + (1 - alpha) * numpy.eye(d) > rs = numpy.random.RandomState(587482) > rv1 = rs.multivariate_normal(mu, R, 1) > print rv1, rv1.sum() > > numpy.random.seed(587482) > rv2 = numpy.random.multivariate_normal(mu, R, 1) > print rv2[0][0], rv2.sum() > > > running it a few times > >>>> > [[ 0.43028555 1.06584226 -0.03496681 -0.31372591 -0.49716804 -1.50641838 > 0.99209124 0.57236839 -0.32107663 0.5865379 ]] 0.973769580498 > 0.0979690286727 0.973769580498 >>>> > [[ 0.09796903 1.41010513 -1.10250773 0.71321445 0.09903517 0.36432555 > -1.27590062 0.04533834 0.37426153 0.24792873]] 0.973769580498 > 0.430285553017 0.973769580498 >>>> > [[ 0.09796903 1.41010513 -1.10250773 0.71321445 0.09903517 0.36432555 > -1.27590062 0.04533834 0.37426153 0.24792873]] 0.973769580498 > 0.430285553017 0.973769580498 >>>> > [[ 0.43028555 1.06584226 -0.03496681 -0.31372591 -0.49716804 -1.50641838 > 0.99209124 0.57236839 -0.32107663 0.5865379 ]] 0.973769580498 > 0.0979690286727 0.973769580498 > > Josef
numpy linalg.svd doesn't produce always the same results running this gives two different answers, using scipy.linalg.svd I always get the same answer, which is one of the numpy answers (numpy random.multivariate_normal is collateral damage) What I don't understand is that numpy.random uses numpy.dual.svd which I thought is scipy.linalg if available, but it looks like it takes the numpy svd. -------- import numpy as np #from numpy.dual import svd from numpy.linalg import svd #from scipy.linalg import svd d = 10 alpha = 1 / d**0.5 mu = numpy.ones(d) R = alpha * numpy.ones((d, d)) + (1 - alpha) * numpy.eye(d) for i in range(10): (u,s,v) = svd(R) print 'v[-1]', v[-1] ----------- Josef _______________________________________________ NumPy-Discussion mailing list NumPy-Discussion@scipy.org http://mail.scipy.org/mailman/listinfo/numpy-discussion