library(ggplot2)
dev.new()
set.seed(1684)
x = seq(15, 33, by = 0.9)
f <- function(x) {
out <- ifelse(
x < 15 | 33 < x,
0,
ifelse(
15 <= x & x <= 24,
(2*(x-15))/((33-15)*(24-15)),
ifelse(
24 < x & x <= 33,
(2*(33-x))/((33-15)*(33-24)),
NA_real_
)))
if (any((is.na(out) | is.nan(out)) & (!is.na(x) & !is.nan(x)))) {
warning("f(x) undefined for some input values")
}
out
}
x_accept <- numeric(0)
while (length(x_accept) != 105) {
x1 = runif(1, min = 15, max = 33)
num = runif(1, min = 0, max = 1)
if (num < (f(x1)/2/(33-15))) {
x_accept = c(x_accept, x1)
}
}
histo <- data.frame(x = x_accept)
dat <- data.frame(x = x, y = f(x), z = histo)
ggplot(dat) +
geom_line(aes(x, y), color = "red") +
geom_histogram(aes(histo), alpha = .5, binwidth = 1)
当我使用这段代码时,它给了我错误,我认为它假设'histo'是一个列表,但它是一个数据框架,所以我真的不知道什么是错的。
使用data=histo
,设置x=x
和y为density (y=..density..
)
ggplot(dat) +
geom_line(aes(x, y), color = "red") +
geom_histogram(aes(x,..density..),data=histo, alpha = .5, binwidth = 1)
或只使用一个数据帧:
dat <- data.frame(x = x, y = f(x), z = x_accept) #z=x_accept instead of histo
ggplot(dat) +
geom_line(aes(x, y), color = "red") +
geom_histogram(aes(z,..density..), alpha = .5, binwidth = 1)