StatRed: added comments to q22.

parent 96dd3447
.PHONY: all clean
all:
clean:
rm -vf *.pyc q*.pdf
from q21_multivariate import dataset from q21_multivariate import dataset
from numpy import array, mean, tile, newaxis, dot from pylab import array, mean, tile, newaxis, dot, eigvals, \
from pylab import eigvals, axis, figure, clf, show, plot axis, figure, clf, show, plot
def eigenvalues(n): def eigenvalues(n):
"""Return eigenvalues of unbiased estimators for the covariance matrix
Sigma (based on a pseudo-random generated dataset)."""
Y = array([mean(dataset(), 1) for i in range(n)]).T Y = array([mean(dataset(), 1) for i in range(n)]).T
# Sigma = 1 / (n - 1) * Sum for i=1 to n: (x_i - x_mean) T(x_i - x_mean),
# where T(x) is the transpose of `x'. Mu = x_mean and
# Yzm = Sum for i=1 to n: x_i - x_mean.
mu = mean(Y, 1) mu = mean(Y, 1)
Yzm = Y - tile(mu[:,newaxis], n) Yzm = Y - tile(mu[:,newaxis], n)
S = dot(Yzm, Yzm.T) / (n - 1) S = dot(Yzm, Yzm.T) / (n - 1)
...@@ -11,14 +18,17 @@ def eigenvalues(n): ...@@ -11,14 +18,17 @@ def eigenvalues(n):
figure(1) figure(1)
clf() clf()
samples = range(2, 10000, 500)
max_range = 10000
samples = range(2, max_range, 500)
data = [[] for i in range(4)] data = [[] for i in range(4)]
for n in samples: for n in samples:
e = eigenvalues(n) e = eigenvalues(n)
for i in range(4): for i in range(4):
data[i].append(e[i]) data[i].append(e[i])
for i in range(4): for i in range(4):
#subplot(2, 2, i+1)
plot(samples, data[i], 'x') plot(samples, data[i], 'x')
axis([0, 10000, 0., 0.025]) axis([0, max_range, 0., 0.025])
show() show()
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