Moysey Abramowitz
Moysey Abramowitz

Reputation: 380

Vectorizing R code, which already uses sequence notation

I am struggling to write fast code to compute a function of the following vector:

enter image description here

Currently I code it using for loop, which is very slow:

rho <- 0.9
E_D <- numeric(100)
E_D[1] <- 1
for (t in 2:100){
  summm <- sum(cumsum(0.9^(0:(t-2)))^2)
  E_D[t] <- t+exp(summm)
}

summm is the element of the vector I analytically defined in the picture above. E_D is a vector, which is some function of that vector. If I set maximum t to 5000, then the code above runs for more than 1 sec on my machine, which is too slow for my purposes.

I tried data.table solution, but it can not accommodate intermediate vector output within a cell:

tempdt <- data.table(prd=2:100 ,summm=0)
tempdt[, summm:=sum(cumsum(rho^(0:(prd-2)))^2)]
Warning message:
In 0:(prd - 2) : numerical expression has 99 elements: only the first used

How to make the code above faster? Please do not tell me that I have to do it in Matlab...

EDIT: To clarify, I need to compute the following vector:

enter image description here

Upvotes: 2

Views: 98

Answers (2)

chinsoon12
chinsoon12

Reputation: 25225

Maybe something like:

n <- 100L
cp <- cumprod(rep(0.9, n - 1L)) / 0.9
cssq <- cumsum(cp)^2
cumsum(cssq)

truncated output:

[1]    1.00000    4.61000   11.95410   23.78082   40.55067   62.50542   89.72283  122.15959 ...

Upvotes: 4

akrun
akrun

Reputation: 886948

The : is not vectorized. We may need to either loop with sapply/lapply or do a group by row

library(data.table)
tempdt[, summm := sum(cumsum(rho^(0:(prd-2)))^2), seq_len(nrow(tempdt))]
head(tempdt)
#   prd    summm
#1:   2  1.00000
#2:   3  4.61000
#3:   4 11.95410
#4:   5 23.78082
#5:   6 40.55067
#6:   7 62.50542

Upvotes: 1

Related Questions