Does R have a way to say "do n times" without having to write a for loop with a redundant variable?

Viewed 504

I recently wrote this line of code: for(i in seq_len(exponent)){out<-squareMat%*%out}. Clearly, i is never used and I just wanted to say "exponent times, do out<-squareMat%*%out". Was there a way to do this without the redundant counter i? For example, is there an apply family function for this?

Example - I have:

squareMat<-matrix(c(0,10,20,30),2,2)
out<-diag(nrow = nrow(squareMat))
exponent<-5
for(i in seq_len(exponent)){out<-squareMat%*%out}
out

What I want is:

squareMat<-matrix(c(0,10,20,30),2,2)
out<-diag(nrow = nrow(squareMat))
exponent<-5
[do exponent times]{out<-squareMat%*%out}
out
6 Answers

Yes, it does. There exist several useful functions in base R that are not used that often these days. One of them does exactly what you want. The replicate function replicates an expression (expr) n times. It works as follows, let's say that I want to generate 3 different samples of size 5 from a uniform (0 to 1) distribution. It can be easily done using replicate. Take a look at the piece of code below

replicate(n = 3, expr = {runif(n = 5)})
#            [,1]      [,2]       [,3]
# [1,] 0.1944426 0.5158065 0.39892501
# [2,] 0.5676580 0.9940599 0.97385575
# [3,] 0.5570141 0.2274214 0.60239883
# [4,] 0.5074303 0.3526040 0.95445298
# [5,] 0.1931812 0.4593620 0.03283596

The results are automatically organized in an array (matrix in this case). However, you can set the parameter simplify = FALSE. Then, the return will be a list

replicate(n = 3, expr = {runif(n = 5)}, simplify = FALSE)
# [[1]]
# [1] 0.4694347 0.9559887 0.8110113 0.7528089 0.6639614
# 
# [[2]]
# [1] 0.8731027 0.7295846 0.3773571 0.5394776 0.6792322
# 
# [[3]]
# [1] 0.3463870 0.3776352 0.3895620 0.2166284 0.5065204

It is important to notice that each of these replications are independent of each other. If you want something to be replicated sequentially, you have to use either a for loop or another suitable function. There exists, for example, a function called rapply (recursive lapply). However, it has never been clear to me the best way to use it.

For a binary operation like %*% you can use Reduce, optionaly with the init argument (actually here not needed because init is identity matrix, see @nicola's comment, but makes answer more general) :

squareMat<-matrix(c(0,10,20,30),2,2)
exponent<-5
Reduce('%*%', init = diag(nrow = nrow(squareMat)), replicate(exponent, squareMat, simplify=F))
#>          [,1]     [,2]
#> [1,]  7800000 27800000
#> [2,] 13900000 49500000

This can be sped up avoiding the call to replicate by not using the second argument in the calculation (1:exponent only acts as a counter) :

Reduce(function(x,notused) {squareMat %*% x}, init = diag(nrow = nrow(squareMat)), 1:exponent)
#>          [,1]     [,2]
#> [1,]  7800000 27800000
#> [2,] 13900000 49500000

Note that Reduce stays in this case less efficient than a simple loop, so, interesting but probably not optimal:

microbenchmark::microbenchmark(test_reduce(),test_reduce2(),test_loop())
Unit: microseconds
           expr  min   lq   mean median    uq  max neval cld
  test_reduce() 17.2 17.7 21.461   18.0 18.60 97.4   100   c
 test_reduce2()  7.9  8.4  9.463    8.6  8.85 25.9   100  b 
    test_loop()  3.0  3.2  3.786    3.3  3.40 26.2   100 a  

An alternative would be to use Recall() within a function call.

my_mat_power = function(n) {
  if (n == 1) 
    squareMat 
  else 
    Recall(n - 1) %*% squareMat
}

my_mat_power(5)

##          [,1]     [,2]
## [1,]  7800000 27800000
## [2,] 13900000 49500000

If you want to make your matrix a variable, you can make small modifications:

my_mat_power = function(n, mat) {
  if (n == 1) 
    mat 
  else 
    Recall(n - 1, mat) %*% mat
}

Here's an example using sapply.

myfunc <- function(s) {
    return('hello')
}
    
sapply(seq_len(5), myfunc)
    
# [1] "hello" "hello" "hello" "hello" "hello"

let's have a look at this :

library(expm)
squareMat %^% 5

R has almost always a way to avaoid loops. In your case:

library(expm)
out <- (squareMat%^%exponent)%*%out
Related