Auto regressive model

Optimizing R Code with Rcpp

Romain François

Consulting Datactive, ThinkR

Auto regressive model, AR

ar <- function(n, phi, sd){
    x  <- epsilon <- rnorm(n, sd = sd)
    np <- length(phi)

    for( i in seq(np+1, n)){
        x[i] <- sum(x[seq(i-1, i-np)] * phi) + epsilon[i]
    }
    x
}

AR in C++

First 🔃, to fill the np first values

NumericVector x(n) ;
// initial loop
for( ___ ; __ < np ; ___ ){
    x[i] = R::rnorm(___) ; }

Main part with outer and inner 🔃

// outer loop
for( ___ ; ___ ; ___ ){    
    double value = rnorm(___) ;    
    // inner loop
    for( ___ ; ___ ; ___ ){
        value += ___ ;    }
    x[i] = value ;
}

Moving average simulation

ma <- function(n, theta, sd){
    epsilon <- rnorm(n, sd = sd)
    x <- numeric(n)
    nq <- length(theta)

    for( i in seq(nq+1, n)){
        x[i] <- sum(epsilon[seq(i-1, i-nq)] * theta) + epsilon[i]
    }
    x
}

Moving average

#include <Rcpp.h>
using namespace Rcpp ;
// [[Rcpp::export]]
NumericVector ma( int n, double mu, NumericVector theta, double sd ){
    int nq = theta.size() ;

    // generate the noise vector at once
    // using the Rcpp::rnorm function, similar to the R function
    NumericVector eps = Rcpp::rnorm(n, 0.0, sd) ;

    // init the output vector of size n with all 0.0 
    NumericVector x(___) ;

    // start filling the values at index nq + 1
    for( int i=nq+1; i<n; i++){
        ____   }
    return x ;
}

ARMA(p,q) = AR(p) + MA(q)

Let's practice!

Optimizing R Code with Rcpp