RcppArmadillo cheatsheet

Miheer Dewaskar edited this page Jun 14, 2017 · 5 revisions

Hi there

I have been using RcppArmadillo more and more frequently, so thought I would make a cheatsheet/cookbook type reference that translates common R operations into equivalent arma code.

You will need the common header to compile and run the examples. Just copy and paste it into R before running any of the examples. This only needs to be done when you first start.

The functions are all pretty basic and not particularly robust. In particular they do not do any bounds or sanity checking.

You might also enjoy the arma documentation, in particular the matlab/octave syntax conversion example.

There is also an excellent book Seamless R and C++ Integration with Rcpp

Any corrections or additions are most welcome.

Contents

Common Header
rbind & cbind type operations
Subsetting & Selection
Some common operations

Common Header

Run this before any of the examples

library('Rcpp')
library('inline')

rcpp_inc <- '
using namespace Rcpp;
using namespace arma;
'

m1 <- matrix(1:16, nr=4)
m2 <- matrix(17:32, nr=4)
v1 <- 1:10
v2 <- 11:20

rbind & cbind type operations

rbind matrix
cbind matrix
c() vectors
cbind vectors to a matrix
rbind vectors to a matrix

rbind w/matrix

src <- '
mat m1 = as<mat>(m1in);
mat m2 = as<mat>(m2in);
mat out = join_cols(m1, m2);
return(wrap(out));
'

fn <- cxxfunction(signature(m1in="numeric", m2in="numeric"), src, plugin='RcppArmadillo', rcpp_inc)
res <- fn(m1, m2)
test <- rbind(m1, m2)
all.equal(test, res)

cbind w/matrix

src <- '
mat m1 = as<mat>(m1in);
mat m2 = as<mat>(m2in);
mat out = join_rows(m1, m2);
return(wrap(out));
'

fn <- cxxfunction(signature(m1in="numeric", m2in="numeric"), src, plugin='RcppArmadillo', rcpp_inc)
res <- fn(m1, m2)
test <- cbind(m1, m2)
all.equal(test, res)

join two vectors, c()

src <- '
vec v1 = as<vec>(v1in);
vec v2 = as<vec>(v2in);
vec out = join_cols(v1, v2);
//returns a nx1 matrix
return(wrap(out));
'
fn <- cxxfunction(signature(v1in="numeric", v2in="numeric"), src, plugin='RcppArmadillo', rcpp_inc)

res <- fn(v1, v2)
test <- c(v1, v2)
all.equal(test, as.vector(res))

cbind vectors to a matrix

src <- '
vec v1 = as<vec>(v1in);
vec v2 = as<vec>(v2in);
mat m1;
m1.insert_cols(m1.n_cols, v1);
m1.insert_cols(m1.n_cols, v2);
//return(wrap(m1));

//or know size in advance
mat m2(v1.n_elem, 2);
m2.col(0) = v1;
m2.col(1) = v2;

return(wrap(m2));
'
fn <- cxxfunction(signature(v1in="numeric", v2in="numeric"), src, plugin='RcppArmadillo', rcpp_inc)

res <- fn(v1, v2)
test <- cbind(v1, v2)
all.equal(test, res, check.attributes=F)

rbind vectors to a matrix

src <- '
vec v1 = as<vec>(v1in);
vec v2 = as<vec>(v2in);

mat m1;
//vecs are column vectors so transpose
m1.insert_rows(m1.n_rows, v1.t());
m1.insert_rows(m1.n_rows, v2.t());

//return(wrap(m1));

mat m2(2, v1.n_elem);
m2.row(0) = v1.t();
m2.row(1) = v2.t();
return(wrap(m2));
'
fn <- cxxfunction(signature(v1in="numeric", v2in="numeric"), src, plugin='RcppArmadillo', rcpp_inc)

res <- fn(v1, v2)
test <- rbind(v1, v2)
all.equal(test, res, check.attributes=F)

Subsetting & Selection

Arma indexes start from 0, subtract 1 from indexes passed from R to arma, and add 1 to indexes returned from arma to R.

Subset vector, contiguous
Subset vector, non-contiguous
Subset matrix by rows, contiguous
Subset matrix by rows, non-contiguous
Subset matrix by columns, contiguous
Subset matrix by columns, non-contiguous
Subset matrix by row & column, contiguous
Subset matrix by row & column, non-contiguous

Subset vector, contiguous

src <- '
vec v1 = as<vec>(v1in);
int from = as<int>(fromin) - 1;
int to = as<int>(toin) - 1;
vec subset = v1.subvec(from, to);
return(wrap(subset));
'
fn <- cxxfunction(signature(v1in="numeric", fromin="numeric", toin="numeric"), src, plugin='RcppArmadillo', rcpp_inc)

res <- fn(v1, 3, 7)
test <- v1[3:7]
all.equal(test, as.vector(res))

Subset vector, non-contiguous

src <- '
vec v1 = as<vec>(v1in);
uvec idx = as<uvec>(idxin) - 1;
vec subset = v1.elem(idx);
return(wrap(subset));
'
fn <- cxxfunction(signature(v1in="numeric", idxin="numeric"), src, plugin='RcppArmadillo', rcpp_inc)

idx <- sample.int(length(v1), 5)
res <- fn(v1, idx)
test <- v1[idx]
all.equal(test, as.vector(res))

Subset matrix by rows, contiguous

src <- '
mat m1 = as<mat>(m1in);
int from = as<int>(fromin) - 1;
int to = as<int>(toin) - 1;
mat s1 = m1.rows(from, to);
//return(wrap(s1));

mat s2 = m1(span(from, to), span::all);
return(wrap(s2));
'
fn <- cxxfunction(signature(m1in="numeric", fromin="numeric", toin="numeric"), src, plugin='RcppArmadillo', rcpp_inc)

res <- fn(m1, 2, 3)
test <- m1[2:3,]
all.equal(test, res)

Subset matrix by rows, non-contiguous

src <- '
mat m1 = as<mat>(m1in);
uvec idx = as<uvec>(idxin) - 1;
mat s1 = m1.rows(idx);
return(wrap(s1));
'
fn <- cxxfunction(signature(m1in="numeric", idxin="numeric"), src, plugin='RcppArmadillo', rcpp_inc)

m4 <- matrix(1:30, nc=3)
idx <- sample.int(nrow(m4), 5)
res <- fn(m4, idx)
test <- m4[idx,]
all.equal(test, res)

Subset matrix by columns, contiguous

src <- '
mat m1 = as<mat>(m1in);
int from = as<int>(fromin) - 1;
int to = as<int>(toin) - 1;
mat s1 = m1.cols(from, to);
//return(wrap(s1));

mat s2 = m1(span::all, span(from, to));
return(wrap(s2));
'
fn <- cxxfunction(signature(m1in="numeric", fromin="numeric", toin="numeric"), src, plugin='RcppArmadillo', rcpp_inc)

res <- fn(m1, 2, 3)
test <- m1[,2:3]
all.equal(test, res)

Subset matrix by columns, non-contiguous

src <- '
mat m1 = as<mat>(m1in);
uvec idx = as<uvec>(idxin) - 1;
mat s1 = m1.cols(idx);
return(wrap(s1));
'
fn <- cxxfunction(signature(m1in="numeric", idxin="numeric"), src, plugin='RcppArmadillo', rcpp_inc)

m4 <- matrix(1:30, nr=3)
idx <- sample.int(ncol(m4), 5)
res <- fn(m4, idx)
test <- m4[,idx]
all.equal(test, res)

Subset matrix by row & column, contiguous

src <- '
mat m1 = as<mat>(m1in);
int row_from = as<int>(row_fromin) - 1;
int row_to = as<int>(row_toin) - 1;
int col_from = as<int>(col_fromin) - 1;
int col_to = as<int>(col_toin) - 1;

mat s1 = m1.submat(row_from, col_from, row_to, col_to);
//return(wrap(s1));

mat s2 = m1(span(row_from, row_to), span(col_from, col_to));
return(wrap(s2));
'
fn <- cxxfunction(signature(m1in="numeric", row_fromin="numeric", row_toin="numeric",
col_fromin="numeric", col_toin="numeric"), src, plugin='RcppArmadillo', rcpp_inc)

m4 <- matrix(1:100, nr=10)
res <- fn(m4, 2, 5, 6, 9)
test <- m4[2:5, 6:9]
all.equal(test, res)

Subset matrix by row & column, non-contiguous

src <- '
mat m1 = as<mat>(m1in);
uvec rowidx = as<uvec>(rowidx_in) - 1;
uvec colidx = as<uvec>(colidx_in) - 1;

mat s1 = m1.submat(rowidx, colidx);
return(wrap(s1));
'
fn <- cxxfunction(signature(m1in="numeric", rowidx_in="numeric", colidx_in="numeric"), src, plugin='RcppArmadillo', rcpp_inc)

m4 <- matrix(1:100, nr=10)
rows = sample.int(nrow(m4), 3);
cols = sample.int(ncol(m4), 3);
res <- fn(m4, rows, cols)
test <- m4[rows, cols]
all.equal(test, res)

Some common operations

rep
diag
bdiag
inverse
crossprod one matrix
crossprod two matrix
tcrossprod one matrix
tcrossprod two matrix
determinant
as.vector, matrix to vector
vector to matrix

rep

src <- '
double what = as<double>(whatin);
int n = as<int>(nin);
vec v(n);
v.fill(what);
return(wrap(v));
'
fn <- cxxfunction(signature(whatin="numeric", nin="numeric"), src, plugin='RcppArmadillo', rcpp_inc)

res <- fn(pi, 10)
test <- rep(pi, 10)
all.equal(test, as.vector(res))

diag

src <- '
vec v1 = as<vec>(v1in);

mat m1 = diagmat(v1);
return(wrap(m1));
'
fn <- cxxfunction(signature(v1in="numeric"), src, plugin='RcppArmadillo', rcpp_inc)

res <- fn(v1)
test <- diag(v1)
all.equal(test, res)

block diag/bdiag

src <- '
vec v1 = as<vec>(v1in);
vec v2 = as<vec>(v2in);

mat m1 = diagmat(join_cols(v1, v2));
return(wrap(m1));
'
fn <- cxxfunction(signature(v1in="numeric", v2in="numeric"), src, plugin='RcppArmadillo', rcpp_inc)

res <- fn(v1, v2)
test <- diag(c(v1, v2))
all.equal(test, res)

inverse

There is also pinv() for the pseudo-inverse

src <- '
mat m1 = as<mat>(m1in);
//return(wrap(m1.i()));

mat m2 = inv(m1);
return(wrap(m2));
'
fn <- cxxfunction(signature(m1in="numeric"), src, plugin='RcppArmadillo', rcpp_inc)

m3 <- matrix(runif(16), nr=4)
res <- fn(m3)
test <- solve(m3)
all.equal(test, res)

crossprod one matrix

src <- '
mat m1 = as<mat>(m1in);
mat cp = trans(m1) * m1;
return(wrap(cp));
'
fn <- cxxfunction(signature(m1in="numeric"), src, plugin='RcppArmadillo', rcpp_inc)

res <- fn(m1)
test <- crossprod(m1)
all.equal(test, res)

crossprod two matrix

src <- '
mat m1 = as<mat>(m1in);
mat m2 = as<mat>(m2in);
mat cp = trans(m1) * m2;
return(wrap(cp));
'
fn <- cxxfunction(signature(m1in="numeric", m2in="numeric"), src, plugin='RcppArmadillo', rcpp_inc)

res <- fn(m1, m2)
test <- crossprod(m1, m2)
all.equal(test, res)

tcrossprod one matrix

src <- '
mat m1 = as<mat>(m1in);
mat tcp = m1 * trans(m1);
return(wrap(tcp));
'
fn <- cxxfunction(signature(m1in="numeric"), src, plugin='RcppArmadillo', rcpp_inc)

res <- fn(m1)
test <- tcrossprod(m1)
all.equal(test, res)

tcrossprod two matrix

src <- '
mat m1 = as<mat>(m1in);
mat m2 = as<mat>(m2in);
mat tcp = m1 * trans(m2);
return(wrap(tcp));
'
fn <- cxxfunction(signature(m1in="numeric", m2in="numeric"), src, plugin='RcppArmadillo', rcpp_inc)

res <- fn(m1, m2)
test <- tcrossprod(m1, m2)
all.equal(test, res)

determinant

src <- '
mat m1 = as<mat>(m1in);
double val, sign;
log_det(val, sign, m1);
return(List::create(Named("modulus")=val, Named("sign") = sign));
'
fn <- cxxfunction(signature(m1in="numeric"), src, plugin='RcppArmadillo', rcpp_inc)

m3 <- matrix(runif(16), nr=4)
res <- unlist(fn(m3))
test <- unlist(determinant(m3))
all.equal(test, res)

as.vector, matrix to vector

src <- '
mat m1 = as<mat>(m1in);
vec v = vectorise(m1);
return(wrap(v));
'
fn <- cxxfunction(signature(m1in="numeric"), src, plugin='RcppArmadillo', rcpp_inc)

res <- fn(m1) #n x 1 matrix
test <- as.vector(m1)
all.equal(test, as.vector(res))

vector to matrix

src <- '
vec v1 = as<vec>(v1in);
int nrow = as<int>(nrow_in);
int ncol = as<int>(ncol_in);
mat m1;
m1.insert_cols(0, v1);
m1.reshape(nrow, ncol);
return(wrap(m1));
'
fn <- cxxfunction(signature(v1in="numeric", nrow_in="numeric", ncol_in="numeric"), src, plugin='RcppArmadillo', rcpp_inc)

res <- fn(1:100, 10, 10)
test <- matrix(1:100, nr=10, nc=10)
all.equal(test, res)

That's it for now!

Clone this wiki locally
You can’t perform that action at this time.
You signed in with another tab or window. Reload to refresh your session. You signed out in another tab or window. Reload to refresh your session.
Press h to open a hovercard with more details.