function [X,acc] = MHsimple
% MHSimple - Simple Metropolis-Hastings Example
% On input:
%     N/A
% On output:
%     X (nx1 vector): samples
%     acc (nx2 array): record of sample acceptance
% Call:
%     [X,acc] = MHSimple;
% Author:
%     T. Henderson
%     UU
%     Winter 2019
%     after the work of Navarro and Perfors
% www.psychology.adelaide.edu.au/personalpages/staff/danielnavarro/ccs2011.html

burnin = 0;
lag = 1;
nsamp = 1000;
sig = 1;
x = -1;

X = zeros(nsamp,1);
acc = [0 0];
for i = 1:burnin
    [x,a] = MHstep(x,sig);
    acc = acc + [a 1];
end
for i = 1:nsamp
    for j = 1:lag
        [x,a] = MHstep(x,sig);
        acc = acc + [a 1];
    end
    X(i) = x;
end

function [x1,a] = MHstep(x0,sig)
%

xp = normrnd(x0,sig);
accprob = targetdist(xp)/targetdist(x0);
u = rand;
if u <= accprob
    x1 = xp;
    a = 1;
else
    x1 = x0;
    a = 0;
end

function probX = targetdist(x)

probX = exp(-x.^2).*(2 + sin(x*5) + sin(x*2));
