Commit 337aed09 authored by Taddeüs Kroes's avatar Taddeüs Kroes

Merge branch 'master' of ssh://vo20.nl/home/git/repos/uva

parents 3f0c473f f3c4b14d
#from math import mean
#import matplotlib.pyplot as plt
from numpy.random import multivariate_normal
def mu_est(x):
return sum(x) / float(len(x))
def cov_est(x):
avg = mean(x)
n = len(x)
# TODO: transpose the second "(i - avg)"
return sum([(i - avg)*(i - avg) for i in x]) / float(n-1)
mean = [0,0]
cov = [[1,0],[0,100]] # diagonal covariance, points lie on x or y-axis
N = 512 # Samples to use.
x,y = multivariate_normal(mean, cov, N).T
print 'x:', x
print 'y:', y
print 'mu_est(x):', mu_est(x)
print 'mu_est(y):', mu_est(y)
mean = [0,0]
cov = [[1,0],[0,100]] # diagonal covariance, points lie on x or y-axis
import matplotlib.pyplot as plt
import numpy as np
x,y = np.random.multivariate_normal(mean, cov, 5000).T
print 'x:', x
print 'y:', y
plt.plot(x,y,'x'); plt.axis('equal'); plt.show()
from numpy import pi, exp, power, dot, rank
from numpy.linalg import inv, det
def mvnd(x, mu, S):
"""
Calculate multivariate normal distribution of vector X, given mu (mean of
elements in vector X) and S (covariance matrix).
Based on the following formula:
http://upload.wikimedia.org/math/a/0/a/a0ad9db46854c4a616ce6959095cf21d.png
"""
return exp((2 * pi) ** rank(x)) * power(det(S), -.5) \
* exp(-.5 * dot(dot((x - mu).T, inv(S)), (x - mu)))
#mean = [0,0]
#cov = [[1,0],[0,100]] # diagonal covariance, points lie on x or y-axis
#
#import matplotlib.pyplot as plt
#import numpy as np
#x,y = np.random.multivariate_normal(mean,cov,5000).T
##mvnd(x, mean, cov)
#plt.plot(x,y,'x'); plt.axis('equal'); plt.show()
http://docs.scipy.org/doc/numpy/reference/generated/numpy.random.multivariate_normal.html
http://packages.python.org/algopy/examples/covariance_matrix_computation.html
http://en.wikipedia.org/wiki/Multivariate_normal_distribution
http://lmf-ramblings.blogspot.com/2009/07/multivariate-normal-distribution-in.html
Markdown is supported
0%
or
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment