Mercurial > hg > camir-aes2014
diff toolboxes/FullBNT-1.0.7/KPMstats/mixgauss_prob_test.m @ 0:e9a9cd732c1e tip
first hg version after svn
author | wolffd |
---|---|
date | Tue, 10 Feb 2015 15:05:51 +0000 |
parents | |
children |
line wrap: on
line diff
--- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/toolboxes/FullBNT-1.0.7/KPMstats/mixgauss_prob_test.m Tue Feb 10 15:05:51 2015 +0000 @@ -0,0 +1,111 @@ +function test_eval_pdf_cond_mixgauss() + +%Q = 10; M = 100; d = 20; T = 500; +Q = 2; M = 3; d = 4; T = 5; + +mu = rand(d,Q,M); +data = randn(d,T); +%mixmat = mk_stochastic(rand(Q,M)); +mixmat = mk_stochastic(ones(Q,M)); + +% tied scalar +Sigma = 0.01; + +mu = rand(d,M,Q); +weights = mixmat'; +N = M*ones(1,Q); +tic; [B, B2, D] = parzen(data, mu, Sigma, N, weights); toc +tic; [BC, B2C, DC] = parzenC(data, mu, Sigma, N); toc +approxeq(B,BC) +B2C = reshape(B2C,[M Q T]); +approxeq(B2,B2C) +DC = reshape(DC,[M Q T]); +approxeq(D,DC) + + +return + +tic; [B, B2] = eval_pdf_cond_mixgauss(data, mu, Sigma, mixmat); toc +tic; C = eval_pdf_cond_parzen(data, mu, Sigma); toc +approxeq(B,C) + +return; + + +mu = reshape(mu, [d Q*M]); + +data = mk_unit_norm(data); +mu = mk_unit_norm(mu); +tic; D = 2 -2*(data'*mu); toc % avoid an expensive repmat +tic; D2 = sqdist(data, mu); toc +approxeq(D,D2) + + +% D(t,m) = sq dist between data(:,t) and mu(:,m) +mu = reshape(mu, [d Q*M]); +D = dist2(data', mu'); +%denom = (2*pi)^(d/2)*sqrt(abs(det(C))); +denom = (2*pi*Sigma)^(d/2); % sqrt(det(2*pi*Sigma)) +numer = exp(-0.5/Sigma * D'); +B2 = numer / denom; +B2 = reshape(B2, [Q M T]); + +tic; B = squeeze(sum(B2 .* repmat(mixmat, [1 1 T]), 2)); toc + +tic +A = zeros(Q,T); +for q=1:Q + A(q,:) = mixmat(q,:) * squeeze(B2(q,:,:)); % sum over m +end +toc +assert(approxeq(A,B)) + +tic +A = zeros(Q,T); +for t=1:T + A(:,t) = sum(mixmat .* B2(:,:,t), 2); % sum over m +end +toc +assert(approxeq(A,B)) + + + + +mu = reshape(mu, [d Q M]); +B3 = zeros(Q,M,T); +for j=1:Q + for k=1:M + B3(j,k,:) = gaussian_prob(data, mu(:,j,k), Sigma*eye(d)); + end +end +assert(approxeq(B2, B3)) + +logB4 = -(d/2)*log(2*pi*Sigma) - (1/(2*Sigma))*D; % det(sigma*I) = sigma^d +B4 = reshape(exp(logB4), [Q M T]); +assert(approxeq(B4, B3)) + + + + +% tied cov matrix + +Sigma = rand_psd(d,d); +mu = reshape(mu, [d Q*M]); +D = sqdist(data, mu, inv(Sigma))'; +denom = sqrt(det(2*pi*Sigma)); +numer = exp(-0.5 * D); +B2 = numer / denom; +B2 = reshape(B2, [Q M T]); + +mu = reshape(mu, [d Q M]); +B3 = zeros(Q,M,T); +for j=1:Q + for k=1:M + B3(j,k,:) = gaussian_prob(data, mu(:,j,k), Sigma); + end +end +assert(approxeq(B2, B3)) + +logB4 = -(d/2)*log(2*pi) - 0.5*logdet(Sigma) - 0.5*D; +B4 = reshape(exp(logB4), [Q M T]); +assert(approxeq(B4, B3))