如何从 cdf 不是封闭形式的分布生成数据？

问题描述

我正在开发一个发行版，其 PDF 和 CDF 为

显然，CDF 的形式并不整洁。那么如何从 R 中的这个 CDF 生成数据，因为这里不能应用逆 CDF 方法？ 我正在使用以下 R 程序从这个 dist 生成一个样本，但我无法证明每一步都是合理的。这是我的代码：

alpha=1.5; beta=3.2;

pdf=function(x)((beta/alpha)*((x/alpha)^beta)*sin(pi/beta))/(((1+(x/alpha)^beta)^2)*(pi/beta))

rand_smplfunc <- function(func,lower,upper){
  x_values <- seq(lower,upper,by = 10^-3)
  sample(x = x_values,size = 1,replace = TRUE,prob = func(x_values))
}

replicate(10,rand_smplfunc(pdf,10))

解决方法

您可以通过以下方式进行，使用具有经验 cdf 的逆 cdf 方法（不需要解析导出的 cdf）：

alpha = 1.5
beta = 3.2

pdf = function(x)((beta/alpha)*((x/alpha)^beta)*sin(pi/beta))/(((1+(x/alpha)^beta)^2)*(pi/beta))
cdf = function(x) cumsum(pdf(x)) / sum(pdf(x)) # approximate the cdf

x <- seq(0.1,50,0.1)
y <- pdf(x)

plot(x,pdf(x),type='l',main='pdf',ylab='f(x)') # plot the pdf

plot(x,cdf(x),main='cdf',ylab='F(x)') # plot the cdf

# draw n samples by using probability integral transform from unifrom random 
# samples by computing the inverse cdf approximately
n <- 5000
random_samples <- approx(x = cdf(x),y = x,runif(n))$y

# plot histogram of the samples
hist(random_samples,200)

从上图可以看出，生成的样本分布确实和原始pdf很像。

,

您可以使用经验 cdf 函数 ecdf 作为近似值。例如

set.seed(7)
pdf <- rnorm(10000)
cdf <- ecdf(pdf)
summary(cdf)
Empirical CDF:    10000 unique values with summary
Min.   1st Qu.    Median      Mean   3rd Qu.      Max. 
-3.930902 -0.671235 -0.001501  0.001928  0.669820  3.791824
range <- -3:3
cdf(range)
[1] 0.0013 0.0219 0.1601 0.5007 0.8420 0.9748 0.9988

因此类似地生成足够数量的 pdf 实现，然后使用 ecdf 函数来近似您的 cdf。

请注意，您还可以绘制经验 cdf：plot(cdf)