We have learnt various R profiling techniques in the previous article. In this twentieth article in the R series, we shall learn how to optimise R code.
We will use version 4.2.2 installed on Parabola GNU/Linux-libre (x86-64) for the code snippets.
$ R --version R version 4.2.2 (2022-10-31) -- “Innocent and Trusting” Copyright (C) 2022 The R Foundation for Statistical Computing Platform: x86_64-pc-linux-gnu (64-bit) R is free software and comes with ABSOLUTELY NO WARRANTY. You are welcome to redistribute it under the terms of the GNU General Public License versions 2 or 3. For more information about these matters see https://www.gnu.org/licenses/.
Recursive versus non-recursive functions
Consider a recursive function to compute the nth Fibonacci sequence:
recurFib <- function(n) {
if (n == 0) {0}
else if (n == 1 || n == 2) {1}
else { recurFib(n - 1) + recurFib(n - 2) }
}
You can verify the output sequence numbers, as shown below:
> recurFib(0) [1] 0 > recurFib(1) [1] 1 > recurFib(2) [1] 1 > recurFib(3) [1] 2 > recurFib(4) [1] 3 > recurFib(5) [1] 5
A non-recursive implementation is given by the following fibo() function:
> fibo <- function(n) {
+ phi <- (sqrt(5) + 1)/2
+ fib <- (phi^(n+1) - (1-phi)^(n+1)) / (2*phi - 1)
+ round(fib)
+ }
In this example, we are starting the Fibonacci number with 1, as illustrated below:
> fibo(0) [1] 1 > fibo(1) [1] 1 > fibo(2) [1] 2 > fibo(3) [1] 3 > fibo(4) [1] 5 > fibo(5) [1] 8
The system.time() function can be used to obtain the computation time taken by both these functions.
> recurFib(35) [1] 9227465 > system.time(recurFib(35)) user system elapsed 6.831 0.000 6.834 > fibo(34) [1] 9227465 > system.time(fibo(34)) user system elapsed 0 0 0
The recursive function is definitely slower than the non-recursive implementation. Moreover, R uses faster native implementation for computing power of values (^), and hence its usage is recommended.
Preallocate
You do not need to allocate memory for any object in R. In the following example, the initialisation takes time.
> v1 <- c()
> system.time(for (i in 1:10000000) {v1[i] <- i;})
user system elapsed
2.143 0.120 2.263
The length can be specified for a vector and memory can be pre-allocated, so that initialisation can be sped up, as illustrated below:
> v2 <- c(NA)
> length(v2) <- 10000000
> system.time(for (i in 1:10000000) {v2[i] <- i;})
user system elapsed
0.299 0.020 0.319
Matrix multiplication
Consider two matrices, m1 and m2:
> m1 <-matrix(1:10000, nrow=100) > m2 <-matrix(10001:20000, nrow=100)
You can write your own version of matrix multiplication as follows:
> matrixMul <- function(mat1, mat2) {
rows <- nrow(mat1)
cols <- ncol(mat2)
output <- matrix(nrow = rows, ncol = cols)
for (i in 1:rows) {
for(j in 1:cols) {
vec1 <- mat1[i,]
vec2 <- mat1[,j]
mult1 <- vec1 * vec2
output[i,j] <- sum(mult1)
}
}
return(output)
}
The time taken for multiplying the two matrices, m1 and m2, using the above function is given below:
> system.time(matrixMul(m1,m2)) user system elapsed 0.039 0.000 0.039
The recommended practice is to perform matrix multiplication using the ‘%*%’ operator. The result for the multiplication of matrices m1 and m2 is quite instantaneous.
> system.time(m1*m2) user system elapsed 0 0 0
Bytecode compiler
The cmpfun() function provides an interface to the byte code compiler for R. It takes the following arguments:
| Argument | Description |
| f | A closure |
| options | Compiler options |
The ‘recurFib’ function can be byte-compiled to improve the performance, as shown below:
> library(compiler) > compiled.recurFib <- cmpfun(recurFib) > system.time(compiled.recurFib(34)) user system elapsed 4.296 0.000 4.298
The compiled function, along with the reference to the bytecode, is shown below:
> compiled.recurFib
function(n) {
if (n == 0) {0}
else if (n == 1 || n == 2) {1}
else { recurFib(n - 1) + recurFib(n - 2) }
}
<bytecode: 0x55eff06fcba8>
You can also view the disassembly of the compiled function using the disassemble() function, as follows:
> disassemble(compiled.recurFib)
list(.Code, list(12L, GETVAR.OP, 1L, LDCONST.OP,3L,EQ.OP, 4L,
BRIFNOT.OP, 5L, 13L, LDCONST.OP, 3L,RETURN.OP, GETVAR.OP,
1L, LDCONST.OP, 8L, EQ.OP, 9L, OR1ST.OP,7L,30L, GETVAR.OP,
1L, LDCONST.OP, 10L, EQ.OP, 11L, OR2ND.OP, 7L, BRIFNOT.OP,
12L, 36L, LDCONST.OP, 8L, RETURN.OP, GETFUN.OP, 14L, MAKEPROM.OP,
17L, CALL.OP, 15L, GETFUN.OP, 14L, MAKEPROM.OP, 19L, CALL.OP,
18L, ADD.OP, 20L, RETURN.OP), list({
if (n == 0) {
0
}
else if (n == 1 || n == 2) {
1
}
else {
recurFib(n - 1) + recurFib(n - 2)
}
}, n, structure(c(2L, 3L, 4L, 44L, 3L, 44L, 2L, 4L), srcfile = <environment>, class = “srcref”),
0, n == 0, if (n == 0) {
0
...
The cmpfile() function can be used to compile a large source .R file. It accepts the following arguments:
| Argument | Description |
| infile | A closure |
| outfile | Compiler options |
| ascii | If compiled file should be stored as ascii |
| env | Environment for compiling |
| verbose | Display compilation output |
| options | Compiler options |
| version | Workspace format version |
Bigmemory package
The bigmemory package is used to handle large, Big Data matrices using shared memory and memory-mapped files. You can install the package as follows:
> install.packages(“bigmemory”) Installing package into ‘/home/guest/R/x86_64-pc-linux-gnu-library/4.1’ (as ‘lib’ is unspecified) --- Please select a CRAN mirror for use in this session --- ... Creating a generic function for ‘typeof’ from package ‘base’ in package ‘bigmemory’ ** help *** installing help indices ** building package indices ** testing if installed package can be loaded from temporary location ** checking absolute paths in shared objects and dynamic libraries ** testing if installed package can be loaded from final location ** testing if installed package keeps a record of temporary installation path * DONE (bigmemory)
The ‘m’ matrix is created using the ‘matrix’ function and written to a .csv file using the ‘write.table()’ function, as shown below:
> library(bigmemory) > m <- matrix(rnorm(30*15000, mean=0, sd=4), 30, 15000) > write.table(m, file = “/tmp/m1.csv”) > system.time(write.table(m, file = “/tmp/m1.csv”)) user system elapsed 0.318 0.001 0.318
You can instead use the big.matrix() function from the bigmemory package to improve the write performance, as demonstrated below:
> system.time(big.matrix(m, nrow=30, ncol=15000, backingfile=”m2.csv”, type=”double”, descriptorfile=”m2.desc”)) user system elapsed 0.006 0.000 0.006
The contents of the m2.csv file are given below for reference:
$ cat m2.desc new(“big.matrix.descriptor”, description = list(sharedType = “FileBacked”, filename = “m2.csv”, dirname = “/home/guest/20-optimization/”, totalRows = 30L, totalCols = 15000L, rowOffset = c(0, 30), colOffset = c(0, 15000), nrow = 30, ncol = 15000, rowNames = NULL, colNames = NULL, type = “double”, separated = FALSE))
The R code optimizer
The rco package takes in R code as input, and applies various optimization techniques to improve the code. You can install the package using the following command:
> install.packages(“rco”) Installing package into ‘/home/guest/R/x86_64-pc-linux-gnu-library/4.1’ (as ‘lib’ is unspecified) --- Please select a CRAN mirror for use in this session --- ... *** installing help indices ** building package indices ** installing vignettes ** testing if installed package can be loaded from temporary location ** testing if installed package can be loaded from final location ** testing if installed package keeps a record of temporary installation path * DONE (rco)
The formula to convert Fahrenheit to Celsius is given below:
celsius = (fahrenheit - 32) / 1.8
The code snippets are stored in a ‘code’ variable, as shown below:
> library(rco)
> code <- paste(
+ “f <- 100”,
+ “paren <- (f - 32)”,
+ “v <- paren / 1.8”,
+ “c <- v”,
+ “if (c > 30) {“,
+ ‘ print(“It is hot!”)’,
+ “} else {“,
+ ‘ print(“It is pleasant!”)’,
+ “}”,
+ sep=”\n”
+ )
You can invoke the optimize_code() function on the code and view the optimised code, as illustrated below:
> pass1 <- optimize_text(code, iterations = 1)
Optimization number 1
> pass1
[1] “f <- 100\nparen <- (100 - 32)\nv <- paren / 1.8\nc <- v\nif (c > 30) {\n print(\”It is hot!\”)\n} else {\n print(\”It is pleasant!\”)\n}”
> pass2 <- optimize_text(pass1, iterations = 1)
Optimization number 1
> pass2
[1] “f <- 100\nparen <- 68\nv <- 68 / 1.8\nc <- v\nif (c > 30) {\n print(\”It is hot!\”)\n} else {\n print(\”It is pleasant!\”)\n}”
The ‘optimize_text’ function accepts the following arguments:
| Argument | Description |
| text | Code in character vectors |
| optimizers | List of optimizer functions |
| iterations | Number of times the optimizer should run |
You are encouraged to read the ‘Help’ pages of the above R functions to learn more about their arguments, options and usage.
