 58
41  54
52  60
60  75
34  47
28  48
35  33];

params.mean0 = [50 50]';
params.covMat0 = [625 312.5; 312.5 625];
params.nu0 = 4;
params.S0 = params.covMat0;

setSeed(0);
nSamples = 1000;
<span class="comment">% data should be D*N</span>
samples = gibbsSampler(data', params, nSamples);
burnin = 250;
samples.mu = samples.mu(:,burnin:end);
samples.Sigma = samples.Sigma(:,:,burnin:end);
Nsamples = size(samples.mu, 2);

<span class="comment">% Plot posterior on mu</span>
[p, grid1, grid2] = ksdensity2d(samples