我正在使用拒绝方法模拟数据,其中的密度函数由X
给出f(x)= C * e^(x) for all x in [0,1]
。我定义g(x) = 1
和C
是其中的最大值f(x)
等于1/(e-1)
。
我使用以下代码来模拟数据:
rejection <- function(f, C, g, rg, n) {
naccepts <- 0
result.sample <- rep(NA, n)
while (naccepts < n) {
y <- rg(1)
u <- runif(1)
if ( u <= f(y) / (C*g(y)) ) {
naccepts <- naccepts + 1
result.sample[naccepts] = y
}
}
result.sample
}
f <- function(x) ifelse(x>=0 & x<=1, (exp(x)), 0)
g <- function(x) 1
rg <- runif
C <- 1/(exp(1) -1)
result <- rejection(f, C, g,rg, 1000)
然后,我使用histogram
将模拟数据与curve
原始数据pdf
进行比较
hist(result,freq = FALSE)
curve(f, 0, 1, add=TRUE)
但是结果的情节有点奇怪!情节就在这里,所以我正在寻找任何帮助来澄清我的工作中出了什么问题。