user3187955
user3187955

Reputation: 13

replacement has length zero in nested loop

I have a problem with a nested loop, here's the code:

T=2
dt=0.001
K=floor(T/dt)
t= seq(0,T,dt)

BM<- matrix(0,2000,1000)
for (i in 1:1000){
for (k in 1:K-1){
BM[k+1,i]=BM[k,i]+sqrt(dt)*rnorm(1)
}
}

the error message is:

Error in BM[k + 1, i] = BM[k, i] + sqrt(dt) * rnorm(1) : 
replacement has length zero

What is wrong in the loop? The goal of it is to simulate 1000 sample paths of Brownian motion.

Upvotes: 1

Views: 2377

Answers (2)

Jilber Urbina
Jilber Urbina

Reputation: 61154

Here's a fully vectorized approach which gives the same result and it is faster than using for loops

B <-  mat.or.vec(K,1000)  # pre-allocate 
set.seed(1)
B[-1, ] <- matrix(rnorm((K-1)*100),K-1,1000) * sqrt(dt)
B <- apply(B, 2, cumsum)

Vectorized functions speed up your code as you can see in the benchmark results below:

          test elapsed relative
1 vectorized()   2.070    1.000
2   for.loop() 181.768   87.811

The vectorized approach is about 88 times faster than your two nested for-loops. Codes for the benchmark are shown below.

for.loop <- function(){  # this is your approach
  BM<- matrix(0,2000,1000)
  set.seed(1)
  for (i in 1:1000){
    for (k in 1:(K-1)){
      BM[k+1, i]=BM[k,i]+sqrt(dt)*rnorm(1)
    }
  }
  return(BM)
}


vectorized <- function(){  # this is my approach
  B <-  mat.or.vec(K,1000)  # pre-allocate
  set.seed(1)
  B[-1, ] <- matrix(rnorm((K-1)*100),K-1,1000) * sqrt(dt)
  B <- apply(B, 2, cumsum)
}


library("rbenchmark")
benchmark(vectorized(),
          for.loop(),
          replications=10,
          columns=c('test', 'elapsed', 'relative'),
          order = "relative")

Comparing the first 5 rows and 5 cols to check results:

> vectorized()[1:5, 1:5]
            [,1]         [,2]         [,3]         [,4]         [,5]
[1,]  0.00000000  0.000000000  0.000000000  0.000000000  0.000000000
[2,] -0.01981021 -0.009865464 -0.008371525  0.015834886  0.001562361
[3,] -0.01400290 -0.037887974 -0.037253019  0.010349073 -0.024924006
[4,] -0.04042779 -0.098675011 -0.073133176  0.018483702 -0.066867547
[5,]  0.01001941 -0.047455576 -0.048955758 -0.001085329 -0.123356698
> for.loop()[1:5, 1:5]
            [,1]         [,2]         [,3]         [,4]         [,5]
[1,]  0.00000000  0.000000000  0.000000000  0.000000000  0.000000000
[2,] -0.01981021 -0.009865464 -0.008371525  0.015834886  0.001562361
[3,] -0.01400290 -0.037887974 -0.037253019  0.010349073 -0.024924006
[4,] -0.04042779 -0.098675011 -0.073133176  0.018483702 -0.066867547
[5,]  0.01001941 -0.047455576 -0.048955758 -0.001085329 -0.123356698

Upvotes: 3

Matthew Lundberg
Matthew Lundberg

Reputation: 42629

This is the problem:

for (k in 1:K-1){

Your precedence is off. This will work:

for (k in 1:(K-1)){

Equivalently (for K > 0):

for (k in seq(K-1)){

Vectorization is better yet:

for (k in seq(K-1)) {
  BM[k+1,] <- BM[k,]+sqrt(dt)*rnorm(2000)
}

Upvotes: 2

Related Questions