%
%  Gibbs-sampler fyrir N(mu,sigma^2)
%
n=500;
mu=1;
sigma=2;
x=randn(n,1)*sigma+mu;
ss=sum((x-mean(x)).^2);
xbar=mean(x);
%
%  Set upp prior
%
mu0=0;
tau2=100000;
a0=1;
b0=1;


nnn=10000;

muvec=zeros(nnn,1);
vvec=zeros(nnn,1);

v=1;
for ii=1:nnn
  b=(n*tau2*xbar+mu0)/(n*tau2+1);
  c=(tau2/(n*tau2*v+v))^0.5;
  mu=randn(1,1)*c+b;
  muvec(ii)=mu;
    bx=1/(ss/2+n*(xbar-mu)^2/2+1/b0);
    v=bx*randg(n/2+a0,1);
    vvec(ii)=v;
end
%
%
%
mean(1./vvec)
var(x)
