Random Walk MH
We can write the following Julia code which use uniform distribution
U[δ,δ]{\mathcal U}[-\delta,\delta]
as
gg
.
1
# random walk Metropolis-Hastings
2
function rmh(T, delta, f::Function, initval = 5)
3
x = ones(T+1)
4
x[1] = initval
5
for t = 1:T
6
# generate Yt
7
epsi = rand() * 2 * delta - delta
8
Yt = epsi + x[t]
9
# accept or not
10
u = rand()
11
r = f(Yt)/f(x[t])
12
if r >= 1
13
x[t+1] = Yt
14
else
15
if u <= r
16
x[t+1] = Yt
17
else
18
x[t+1] = x[t]
19
end
20
end
21
end
22
return(x)
23
end
Copied!
Then apply this algorithm to the normal distribution:
1
# density function of N(0, 1) without normalization
2
function dnorm(x)
3
return(exp(-0.5 * x^2))
4
end
5
6
# example
7
deltalist = [0.1, 0.5, 1]
8
N = 15000
9
# results of mean and variance
10
muvar = zeros(3, 2)
11
using Statistics
12
ps = []
13
for i = 1:3
14
x = rmh(N, deltalist[i], dnorm)
15
push!(ps, plot(x, legend=:none))
16
muvar[i, 1] = mean(x)
17
muvar[i, 2] = var(x)
18
end
19
# plot
20
plot(ps[1], ps[2], ps[3], layout=(3,1))
21
savefig("res_rmh.png")
Copied!
We will get the table of mean and variance as showed in Table 6.3.2 of Robert and Casella (2013)
and the curve of each case:
Last modified 2yr ago
Copy link