##############
# Perform Metropolis Hastings with T iterations 
# to obtain a sample from a Poisson(44) distribution
# Note the steps included in the for loop involve a 
# random proposed value, then a random decision
#
# Stat 380 Week 6 Day 2
# Dr. Dave Campbell
# Simon Fraser University
# Feb 19, 2015
#############

T=10000           #Use this many iterations
x<-rep(0,1,T)     #this makes a vector of zeros to fill with our MCMC values

x[1]<-44          #we need to initialize the method somewhere, may as well be here


for(lp in 1:T){   #I use ‘lp’ as my index name      

	#propose a value for Y from Q where all available states have the same probability
	range<-(x[lp]-3) : (x[lp]+3) #brackets are important
	proposedy<-sample(range,1)   #sample discrete uniform over the set of values in range 
		

	#Compute alpha so that we can make a decision about accepting the proposed value
	alpha<-dpois(proposedy,44)/dpois(x[lp],44) #since the transition matrix is symmetric, i.e.: Qxy=Qyx
	#Decide to keep the proposedy value or retain the last accepted value with probability min(alpha, 1)
	u<-runif(1,0,1)
	ifelse(u<alpha,x[lp+1]<-y,x[lp+1]<-x[lp]) 
	#Now go back to the start of the loop
}