So far we’ve had a look at rejection sampling and importance sampling. Here we take a quick look at slice sampling, although rather than implementing it ourselves, we will use the built in Matlab slicesample function.
Using our parameter estimation example, we will use slice sampling to estimate the mean and sigma of some samples from a normal distribution. Here is the code.
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
| true_mean = 0; | |
| true_sigma = 1; | |
| % likelihood_func = @(x, mean, sigma) normpdf(x, mean, sigma); | |
| % the above function to calcalate in matrix form, for speed | |
| likelihood_func = @(params)… | |
| prod(normpdf(repmat(x,[1 numel(params(1))]),… | |
| repmat(params(1), [1 numel(x)])',… | |
| repmat(params(2),[1 numel(x)])' ), 1); | |
| %% generate data | |
| global x | |
| N=20; | |
| x = normrnd(true_mean, true_sigma, [N 1]); | |
| %% DO SLICE SAMPLING | |
| initial_samples = [0,1]; | |
| nSamples = 10^4; | |
| trace = slicesample(initial_samples, nSamples,… | |
| 'burnin', 100,… | |
| 'pdf',likelihood_func); | |
| mean_samples = trace(:,1); | |
| sigma_samples = trace(:,2); | |
| subplot(1,2,1) | |
| plot(mean_samples,sigma_samples,'.') | |
| xlabel('mean') | |
| ylabel('sigma') | |
| hold on | |
| plot(true_mean, true_sigma, 'r.','MarkerSize',5^2) | |
| axis square | |
| subplot(2,2,2) | |
| plot(mean_samples) | |
| ylabel('mean') | |
| subplot(2,2,4) | |
| plot(sigma_samples) | |
| ylabel('sigma') |
Which results in this plot of samples estimating the posterior P(params|data)
