我有以下任务要执行：
生成公式描述的10 ^ 9步骤：
X(0)=0 X(t+1)=X(t)+Y(t)其中
Y(t)是独立的随机变量，分布为N(0,1)。计算t指数的百分比为X(t)的百分比为负数。
我尝试了以下代码：
x<-c(0,0)
z<-0
loop<-10^9
for(i in 2:loop) {
x[1]<-x[2]
x[2]<-x[1]+rnorm(1, 0, 1)
if (x[2]<0) {z<-z+1}
}
但是，它很慢。我怎样才能加快速度？
一种解决方案是使用@ G5W提出的矢量化，但将其分解为更小的块以避免任何内存溢出问题。这为您提供了矢量化解决方案的速度，但通过管理块大小，您可以控制进程使用的内存量。
以下将问题分解为1e + 07的块，通过循环100次，总共获得1e + 09。
在第一个块的末尾，记录低于0的时间百分比和结束点。然后将结束点送入下一个块，并记录低于0的时间百分比和新的结束点。
最后，平均100次运行以使总时间量低于零。在while循环中调用cat是为了监视进度并查看进度，这可以被注释掉。
funky <- function(start, length = 1e+07) {
Y <- rnorm(length)
Z <- cumsum(Y)
c(sum(Z<(-start))/length, (tail(Z, 1) + start))
}
starttime <- Sys.time()
resvect <- vector(mode = "numeric", length = 100)
result <- funky(0)
resvect[1] <- result[1]
i <- 2
while (i < 101) {
cat(result, "\n")
result <- funky(result[2])
resvect[i] <- result[1]
i <- i + 1
}
mean(resvect)
# [1] 0.1880392
endtime <- Sys.time()
elapsed <- endtime - starttime
elapsed
# Time difference of 1.207566 mins
通常，对于这样的问题，您可以使用Rcpp包将函数一对一地转换为C ++。这应该会带来相当大的加速。
一，R版：
random_sum <- function(loop = 1000) {
x<-c(0,0)
z<-0
for(i in 2:loop) {
x[1]<-x[2]
x[2]<-x[1]+rnorm(1, 0, 1)
if (x[2]<0) {z<-z+1}
}
z / loop
}
set.seed(123)
random_sum()
# [1] 0.134
现在的C ++版本：
library("Rcpp")
cppFunction("
double random_sum_cpp(unsigned long loop = 1000) {
double x1 = 0;
double x2 = 0;
double z = 0;
for (unsigned long i = 2; i < loop; i++) {
x1 = x2;
x2 = x1 + Rcpp::rnorm(1)[0];
if (x2 < 0) z = z+1;
}
return z/loop;
}")
set.seed(123)
random_sum_cpp()
# [1] 0.134
为了完整起见，我们还要考虑提出的矢量化版本：
random_sum_vector <- function(loop = 1000) {
Y = rnorm(loop)
sum(cumsum(Y)<0)/loop
}
set.seed(123)
random_sum_vector()
# [1] 0.134
我们看到它为相同的随机种子提供了相同的结果，因此它似乎是一个可行的竞争者。
在基准测试中，C ++版本和矢量化版本的表现相似，矢量化版本在C ++版本上略显优势：
> microbenchmark(random_sum(100000),
random_sum_vector(100000),
random_sum_cpp(100000))
Unit: milliseconds
expr min lq mean median uq max neval
random_sum(1e+05) 184.205588 199.859266 209.220232 205.137043 211.026740 274.47615 100
random_sum_vector(1e+05) 6.320690 6.631704 7.273645 6.799093 7.334733 18.48649 100
random_sum_cpp(1e+05) 8.950091 9.362303 10.663295 9.956996 11.079513 21.30898 100
但是，矢量化版本通过内存和will blow up your memory for long loops.来降低速度.C ++版本几乎不使用内存。
对于10 ^ 9步，C ++版本在我的机器上运行大约2分钟（110秒）。我没试过R版。根据较短的基准，它可能需要大约7个小时。
> microbenchmark(random_sum_cpp(10^9), times = 1)
Unit: seconds
expr min lq mean median uq max neval
random_sum_cpp(10^9) 110.2182 110.2182 110.2182 110.2182 110.2182 110.2182 1
这应该快得多，但十亿可能需要一段时间。用较小的长度值来测试它可能是好的 - 比如10 ^ 6。
length = 10^9
Y = rnorm(length)
sum(cumsum(Y)<0)/length
基于@ user3666197的评论我测试了这个，他是对的。此解决方案适用于较小的数字，但一旦步骤数太大，它就会失败。
我根据OP的代码测试了我的“矢量化”版本。当随机游走的长度是10 ^ 8时，我的代码花了大约7秒钟，OP的代码花了131秒（在我的笔记本电脑上）。但是，当我将长度增加到10 ^ 9时（根据原始问题），我的版本导致了大量的磁盘交换，我不得不杀死进程。该解决方案以OP请求的规模失败。
给定随机源在技术上被构造为否则确定性硬件的能力，以满足对所生成的流的可重复性的要求以及由某个伪随机生成器算法“生成的”随机性的所有条件，这样的源随机性不容易从纯粹的[SERIAL]转换成任何形式的“just”-[CONCURRENT]或真正的[PARALLEL]运作方式。
这就是说，PRG步骤是重新定义纯[SERIAL]代码执行的任何尝试的中心点（阻塞）。
这不会改变（非）负的X(t)值的百分比，但只是确定对于给定的PRG-硬件实现，没有更短的方式，而是纯[SERIAL]生成相互（串行）相关值的序列。
展开“慢”循环或准循环（因为值依然是串行依赖） - 矢量处理（R语言实现功能利用但几乎硬件CPU指令集级别技巧 - 所以不是语言游戏改变者，但有点绕过一些明显缓慢的代码执行构造函数）是最可能发生的事情。
使用向量通常会产生比循环更好的性能。这里的问题非常大（即10 ^ 9）是存储器限制。由于您只对负数指数的最终百分比感兴趣，因此以下方法可行（在10 ^ 9步骤中需要几分钟）。
update_state <- function (curr_state, step_size) {
n <- min(curr_state$counter, step_size)
r <- rnorm(min(curr_state$counter, step_size))
total <- curr_state$cum_sum + cumsum(r)
list('counter' = curr_state$counter - n,
'neg_count' = curr_state$neg_count + length(which(total < 0)),
'cum_sum' = curr_state$cum_sum + sum(r))
}
n <- 10^9
curr_state <- list('counter' = n, 'neg_count' = 0, 'cum_sum' = 0)
step_size <- 10^8
while (curr_state$counter > 0) {
curr_state <- update_state(curr_state = curr_state, step_size = step_size)
}
print(curr_state)
print(curr_state$neg_count/ n)