⚡ Pluto.jl ⚡

Why Metropolis–Hastings Works

Example: Rosenbrock density

131 μs

xxxxxxxxxx
 
using Distributions, Plots, Random, LinearAlgebra,LaTeXStrings

13.0 s

f (generic function with 1 method)

xxxxxxxxxx
 
begin
    x = range(-3,stop=3,length=200)
    y = range(-3,stop=10,length=200)
    N = 1000
    Sampleₙ = 5
    f(x,y) = exp(-((1 - x)^2 + 100(y - x^2)^2) / 20)
end

140 ms

xxxxxxxxxx
 
plot(x,y,(x,y)->f(x,y),st=:contour,label=L"f(x,y)",legend=:bottomleft)

5.8 s

DiagNormal(
dim: 2
μ: [0.0, 0.0]
Σ: [0.1 0.0; 0.0 0.1]
)

xxxxxxxxxx
 
q = MvNormal(zeros(2),Diagonal(ones(2))*0.1)

151 ms

xxxxxxxxxx
 
begin
    X = -1:0.01:2
    Y = -1:0.01:2
    Z = [pdf(q,[i,j]) for i in X, j in Y] 
    plot(X,Y,Z,st=:contour)
end

111 ms

Julia Plots Color Scheme

3.7 μs

xxxxxxxxxx
 
begin
    Z₂ = [pdf(q,[i,j]) for i in X, j in Y] 
    plot(X,Y,Z₂,st=:surface,camera=(80,45),color=:amp)
end

1.1 s

$\begin{aligned} (x^{*}, y^{*}) = \underset{x, y}{argmax f_{m a x}} \\ f_{m a x} (x^{*} = 1, y^{*} = 1) = 1 \end{aligned}$

참고 : Wolframalpha Solution

3.7 μs

xxxxxxxxxx
 
plot(x,y,(x,y)->f(x,y),st=:surface,camera=(30,80),alpha=0.1,color=:heat)

80.8 ms

 
begin   
    plot(x,y,(x,y)->f(x,y),st=:contour)
    for sᵢ in 1:Sampleₙ
        samples = Array{Float64,2}(zeros(N,2))
        samples[1,:] = [rand(Uniform(-3,3)),rand(Uniform(-3,10))]
        udist = Uniform()
        for i in 2:N
            curr = samples[i-1,:]
            prop = curr + rand(q)
            alpha = f(prop...) / f(curr...)     
            if rand(udist) < alpha
                curr = prop
            end
            samples[i,:] = curr
        end 
        plot!(samples[:,1],samples[:,2],marker=(:circle,3), label=L"sample_{%$sᵢ}",
                legend=:bottomleft)
    end
    plot!([1.0],[1.0],st=:scatter,marker=(:star5,15),color=:yellow,label=L"(1.0,1.0)")
    #plot!([1,2],[1,-2],linewidth=2,label="")
    quiver!([1],[1],quiver=([2-1],[-2-1]),color=:red)
    annotate!(3, -2.5, text(L"f_{max}=1\;at\;(1,1)", :red, :right, 9))
end

58.1 ms