arrays - R - avoid nested for loops -


i have following function takes 4 vectors. t vector has given length , 3 other vectors (pga, sa5hz , sa1hz) have given (identical not equal t) lenght.

the output matrix length(t) rows , length(pga) columns.

my code below seems perfect example of not do, however, not figure out way optimize using apply function. can help?

designspectrum <- function (t, pga, sa5hz, sa1hz){    ts <- sa1hz / sa5hz    #by convention, if sa5hz null, set ts 0.   ts[is.nan(ts)] <- 0    res <- matrix(na, nrow = length(t), ncol = length(pga))    (i in 1:nrow(res))   {     (j in 1:ncol(res))     {       res[i,j] <- if(t[i] <= 0) {pga[j]}                   else if (t[i] <= 0.2 * ts[j]) {pga[j] + t[i] * (sa5hz[j] - pga[j]) / (0.2 * ts[j])}                   else if (t[i] <= ts[j]) {sa5hz[j]}                   else sa1hz[j] / t[i]       }   }    return(res) } 

instead of doing double for loop , processing each i , j value separately, use outer function process of them in 1 shot. since you're processing multiple i , j values simultaneously, switch vectorized ifelse statement instead of non-vectorized if , else statements:

designspectrum2 <- function (t, pga, sa5hz, sa1hz) {   ts <- sa1hz / sa5hz   ts[is.nan(ts)] <- 0   outer(1:length(t), 1:length(pga), function(i, j) {     ifelse(t[i] <= 0, pga[j],       ifelse(t[i] <= 0.2 * ts[j], pga[j] + t[i] * (sa5hz[j] - pga[j]) / (0.2 * ts[j]),         ifelse(t[i] <= ts[j], sa5hz[j], sa1hz[j] / t[i])))   }) } identical(designspectrum(t, pga, sa5hz, sa1hz), designspectrum2(t, pga, sa5hz, sa1hz)) # [1] true 

data:

t <- -1:3 pga <- 1:3 sa5hz <- 2:4 sa1hz <- 3:5 

you can see efficiency gains testing on rather large vectors (here i'll use output matrix 1 million entries):

# larger vectors set.seed(144) t2 <- runif(1000, -1, 3) pga2 <- runif(1000, -1, 3) sa5hz2 <- runif(1000, -1, 3) sa1hz2 <- runif(1000, -1, 3)  # runtime comparison all.equal(designspectrum(t2, pga2, sa5hz2, sa1hz2), designspectrum2(t2, pga2, sa5hz2, sa1hz2)) # [1] true system.time(designspectrum(t2, pga2, sa5hz2, sa1hz2)) #    user  system elapsed  #   4.038   1.011   5.042  system.time(designspectrum2(t2, pga2, sa5hz2, sa1hz2)) #    user  system elapsed  #   0.517   0.138   0.652  

the approach outer 10x faster.


Comments

Popular posts from this blog

html - Outlook 2010 Anchor (url/address/link) -

javascript - Why does running this loop 9 times take 100x longer than running it 8 times? -

Getting gateway time-out Rails app with Nginx + Puma running on Digital Ocean -