Factorial Memoization in R

Question

I wrote this function to find a factorial of number

fact <- function(n) {
    if (n < 0){
      cat ("Sorry, factorial does not exist for negative numbers", "\n")
    } else if (n == 0){
      cat ("The factorial of 0 is 1", "\n")
    } else {
    results = 1
    for (i in 1:n){
      results = results * i
    }
    cat(paste("The factorial of", n ,"is", results, "\n"))
    }
}

Now I want to implement Memoization in R. I have Basic idea on R and trying to implement using them. But I am not sure is this way forward. Could you please also elaborate this topic as well. Thank in advance. Memoized Factorial

    fact_tbl <- c(0, 1, rep(NA, 100))
    fact_mem <- function(n){
          stopifnot(n > 0)
          if(!is.na(fib_tbl[n])){
           fib_tbl[n]
    } else {
       fact_tbl[n-1] <<- fac_mem(n-1) * n
     }
   }

   print (fact_mem(4))

Answer 1

First of all, if you need an efficient implementation, use R's factorial function. Don't write it yourself. Then, the factorial is a good exercise for understanding recursion:

myfactorial <- function(n) {
  if (n == 1) return(1)
  n * myfactorial(n-1)
}

myfactorial(10)
#[1] 3628800

With this function memoization is only useful, if you intend to use the function repeatedly. You can implement memoization in R using closures. Hadley explains these in his book.

createMemFactorial <- function() {
  res <- 1
  memFactorial <- function(n) {
    if (n == 1) return(1)

    #grow res if necessary
    if (length(res) < n) res <<- `length<-`(res, n)

    #return pre-calculated value
    if (!is.na(res[n])) return(res[n])

    #calculate new values
    res[n] <<- n * factorial(n-1)
    res[n]
  }
  memFactorial
}
memFactorial <- createMemFactorial()

memFactorial(10)
#[1] 3628800

Is it actually faster?

library(microbenchmark)
microbenchmark(factorial(10),
               myfactorial(10), 
               memFactorial(10))
#Unit: nanoseconds
#             expr  min     lq    mean median     uq   max neval cld
#    factorial(10)  235  264.0  348.02  304.5  378.5  2463   100 a  
#  myfactorial(10) 4799 5279.5 6491.94 5629.0 6044.5 15955   100   c
# memFactorial(10)  950 1025.0 1344.51 1134.5 1292.0  7942   100  b 

Note that microbenchmark evaluates the functions (by default) 100 times. Since we have stored the value for n = 10 when testing the memFactorial, we time only the if conditions and the lookup here. As you can also see, R's implementation, which is mostly written in C, is faster.

A better (and easier) example implements Fibonacci numbers. Here the algorithm itself benefits from memoization.

#naive recursive implementation
fib <- function(n)  {
  if(n == 1 || n == 2) return(1)
  fib(n-1) + fib(n-2)
}

#with memoization
fibm <- function(n)  {
  if(n == 1 || n == 2) return(1)

  seq <- integer(n)
  seq[1:2] <- 1

  calc <- function(n) {
    if (seq[n] != 0) return(seq[n])
    seq[n] <<- calc(n-1) + calc(n-2)
    seq[n]
  }

  calc(n)
}

#try it:
fib(20)
#[1] 6765
fibm(20)
#[1] 6765

#Is memoization faster?
microbenchmark(fib(20),
               fibm(20))
#Unit: microseconds
#     expr      min       lq       mean    median        uq       max neval cld
# fib(20)  8005.314 8804.130 9758.75325 9301.6210 9798.8500 46867.182   100   b
#fibm(20)    38.991   44.798   54.12626   53.6725   60.4035    97.089   100  a