s2<-1;t2<-10;mu<-5;//gaussian prior~N(mu,t2),likelihood~N(theta,s2) 
y<-c(9.37,10.18,9.16,11.60,10.33)//the data
theta<-0; delta2<-2;  THETA<-NULL; set.seed(1)//sample new parameter from the previous parameter according to N(theta, delta)
S<-10000; // the number of iteration
for (s in 1:S)
{
  theta.star<-rnorm(1, theta,sqrt(delta2));

	log.r<-(sum(dnorm(y,theta.star,sqrt(s2),log=TRUE))+dnorm(theta.star,mu,sqrt(t2),log=TRUE))-(sum(dnorm(y,theta,sqrt(s2),log=TRUE))+dnorm(theta,mu,sqrt(t2),log=TRUE))

      if(log(runif(1))<log.r){theta<-theta.star}

	THETA<-c(THETA,theta)
}