%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%  MH sampling example
%
% Samples from posterior for mean (variance is fixed in this example)
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

function [mu] = MetropolisFirstExpl(x, trueMu, trueSigma)
%input
% x: data
% trueMu: true value of mu (this is only used in the graphics)
% trueSigma: true value of sigma

step = 1.5;
iterations = 10000;
burnin = 1000;

muInit=5;
%priorMu = 1;

%%Create space for the stored variables
mu = zeros(1,iterations);
mu(1) = muInit;
phi=1/trueSigma^2;   %use precision parameterization (this is arbitrary)

hold on
for(i=2:iterations)
    
    %%update the mean--------------------------------------------------- 
    proposalMu = mu(i-1) + randn*step/2;
    proposed    = sum( - phi*(x-proposalMu).^2/2);
    currentlyAt = sum( - phi*(x-mu(i-1)).^2/2);  
       
    logratio = proposed-currentlyAt
    
    if(log(rand)<logratio)
        mu(i) = proposalMu;
    else
        mu(i) = mu(i-1);
    end
    %%end the mean update------------------------------------------------
    
    
   %just plotting the results periodically 
   stepplot=100; 
   if(mod(i,stepplot)==0)  
       i
     plot(mu(1:i),'-','color',[200/255,200/255,200/255])
     plot(mu(i-stepplot+1:i),'k-')
     plot(trueMu,'r*')
     pause(.01)
   end
end


%plot all the output
figure(5)
clf;
subplot(2,1,1)
plot(mu)

mu = mu((burnin+1):length(mu));
title('trace of \mu')

subplot(2,1,2)
%normalizedHist(mu,.2)
hist(mu,20)
hold on;
plot([trueMu,trueMu],[0,length(mu)/10],'r','linewidth',2)
title('distribution of \mu')


  
    