R语言舍选法产生随机数

Question

用R语言舍选法产生随机数，当0≤x≤1时，F（x）=(x+x^2)/2,画出频率直方图、核密度函数曲线以及真实的密度函数曲线

Answer 1

引自免费微信小程序：皆我百晓生

在R语言中，要使用舍选法（Inverse Transform Method）生成符合特定累积分布函数（CDF）的随机数，首先你需要计算出逆累积分布函数（ICDF），然后对均匀分布的随机数进行映射。

根据给定的累积分布函数 ( F(x) = \frac{x + x^2}{2} )，我们先求其逆函数 ( F^{-1}(y) )。由于这个CDF不是标准形式，我们可能需要数值方法来近似解逆函数。这里我们将使用二分法来找到 ( F^{-1}(y) ) 的近似值。

接下来，我会展示如何在R语言中实现这一过程，并绘制频率直方图、核密度函数曲线及真实的密度函数曲线。

# 导入必要的库
library(ggplot2)
library(Metrics)

# 定义原始累积分布函数
F <- function(x) {
  (x + x^2) / 2
}

# 利用二分法求解逆累积分布函数的近似值
find_F_inv <- function(p, lower = 0, upper = 1, tol = 1e-6) {
  if (p <= 0 | p >= 1) {
    stop("p must be between 0 and 1")
  }
  
  while ((upper - lower) > tol) {
    mid <- (lower + upper) / 2
    if (F(mid) < p) {
      lower <- mid
    } else {
      upper <- mid
    }
  }
  
  return((lower + upper) / 2)
}

# 生成随机数
n_samples <- 10000
random_p <- runif(n_samples)
random_x <- sapply(random_p, find_F_inv)

# 绘制频率直方图
hist_random <- hist(random_x, breaks = "Scott", freq = FALSE, main = "Frequency Histogram",
                    xlab = expression(F[x]), ylab = "Density")

# 计算核密度估计
density_random <- density(random_x)

# 绘制真实的密度函数曲线
true_density <- function(x) {
  1 / sqrt(2 * pi) * exp(-(x - 1)^2 / 2)
}
x_values <- seq(-3, 3, length.out = 1000)
true_density_line <- data.frame(x = x_values, y = true_density(x))

# 组合并绘制所有图形
ggplot() +
  geom_histogram(data = hist_random$data, aes(x = hist_random$mids), 
                 fill = "lightblue", alpha = 0.5, bins = hist_random$breaks) +
  geom_line(data = density_random, aes(x = x, y = ..scaled..), color = "red") +
  geom_line(data = true_density_line, aes(x = x, y = y), linetype = "dashed", color = "black") +
  scale_y_continuous(NULL, breaks = NULL) +
  theme_classic() +
  labs(title = "Random Numbers and Distributions",
       x = expression(F[x]),
       y = "Density")

这段代码首先定义了累积分布函数 F 和逆累积分布函数的近似求解器 find_F_inv。接着，它生成了1万个遵循给定CDF的随机数，并计算了它们的频率直方图和核密度估计。同时，也绘制了真实的密度函数曲线作为参考。最后，将这三种图形合并到一个ggplot2图中进行可视化。