我正在尝试使用 R 从自定义 pdf 生成随机样本。我的 pdf 是：
我生成了统一的样本，然后尝试将其转换为我的自定义分布。我通过找到我的发行版的 cdf 来做到这一点（) 并将其设置为均匀样本 () 并求解.
生成具有上述分布的随机样本，得到一个均匀的样本并解决在
我实现了它，R但没有得到预期的分布。谁能指出我理解中的缺陷？
nsamples <- 1000;
x <- runif(nsamples);
f <- function(x, u) {
return(3/2*(x-x^3/3) - u);
}
z <- c();
for (i in 1:nsamples) {
# find the root within (0,1)
r <- uniroot(f, c(0,1), tol = 0.0001, u = x[i])$root;
z <- c(z, r);
}
看起来您发现您的代码有效，但@Aniko 指出您可以提高其效率。您最大的速度提升可能来自预分配内存，z这样您就不会在循环中增长它。类似的东西z <- rep(NA, nsamples)应该可以解决问题。vapply()您可能会通过使用（指定返回的变量类型）而不是显式循环（在 apply 系列上有一个很好的SO 问题）获得小的速度增益。
> nsamples <- 1E5
> x <- runif(nsamples)
> f <- function(x, u) 1.5 * (x - (x^3) / 3) - u
> z <- c()
>
> # original version
> system.time({
+ for (i in 1:nsamples) {
+ # find the root within (0,1)
+ r <- uniroot(f, c(0,1), tol = 0.0001, u = x[i])$root
+ z <- c(z, r)
+ }
+ })
user system elapsed
49.88 0.00 50.54
>
> # original version with pre-allocation
> z.pre <- rep(NA, nsamples)
> system.time({
+ for (i in 1:nsamples) {
+ # find the root within (0,1)
+ z.pre[i] <- uniroot(f, c(0,1), tol = 0.0001, u = x[i])$root
+ }
+ })
user system elapsed
7.55 0.01 7.78
>
>
>
> # my version with sapply
> my.uniroot <- function(x) uniroot(f, c(0, 1), tol = 0.0001, u = x)$root
> system.time({
+ r <- vapply(x, my.uniroot, numeric(1))
+ })
user system elapsed
6.61 0.02 6.74
>
> # same results
> head(z)
[1] 0.7803198 0.2860108 0.5153724 0.2479611 0.3451658 0.4682738
> head(z.pre)
[1] 0.7803198 0.2860108 0.5153724 0.2479611 0.3451658 0.4682738
> head(r)
[1] 0.7803198 0.2860108 0.5153724 0.2479611 0.3451658 0.4682738
而且您不需要;每行末尾的 （您是 MATLAB 转换器吗？）。