Open henryliangt opened 1 year ago
knitr::opts_chunk$set(echo = TRUE)
set.seed(2)
p <- 0.4
mn <- c(-1,2)
std <- c(0.2, 1.5)
f <- function(x, mu=mn, sd=std) p * dnorm(x, mu[1], sd[1]) + (1-p) * dnorm(x, mu[2], sd[2])
curve(f(x), col="blue", from=-4, to=8, n=300, las=1)
q <- function(x, sigma) rnorm(1, mean=x, sd=sigma)
mc.algorithm <- function(x_0, n, sigma)
{
x <- numeric(n)
x[1] <- x_0
for (i in 1:(n-1))
{
prop.x <- q(x[i], sigma=sigma)
alpha <- min(f(prop.x)/f(x[i]), 1)
if (runif(1) < alpha)
x[i+1] <- prop.x
else
x[i+1] <- x[i]
}
x
}
choices <- expand.grid(x_0 = c(-1,1,5),sigma=c(0.01,0.5,20))
par(mfrow=c(2,2))
x <- seq(from=-4,to=8,length.out = 300)
dresult <- f(x)
for(i in 1:nrow(choices))
{
result <- mc.algorithm(choices$x_0[i],1000,choices$sigma[i])
hist(result[100:1000],n = 50,freq = F, main = paste0('x_0 = ',choices$x_0[i],' sigma = ',choices$sigma[i]))
lines(x, dresult, col='blue')
}
handbook: https://towardsdatascience.com/how-to-use-r-in-google-colab-b6e02d736497
from here: https://colab.research.google.com/#create=true&language=r
URL https://colab.to/r
R.version.string
print(installed.packages())