Using RcppArmadillo to price European Put Options

By Rcpp Gallery

(This article was first published on Rcpp Gallery, and kindly contributed to R-bloggers)

Introduction

In the quest for ever faster code, one generally begins exploring ways to integrate C++
with R using Rcpp. This post provides an example of multiple
implementations of a European Put Option pricer. The implementations are done in pure R,
pure Rcpp using some Rcpp sugar functions,
and then in Rcpp using
RcppArmadillo, which exposes the
incredibly powerful linear algebra library, Armadillo.

Under the Black-Scholes model The value of a European put option has the closed form solution:

where

and

Armed with the formulas, we can create the pricer using just R.

put_option_pricer  function(s, k, r, y, t, sigma) {

    d1  (log(s / k) + (r - y + sigma^2 / 2) * t) / (sigma * sqrt(t))
    d2  d1 - sigma * sqrt(t)

    V  pnorm(-d2) * k * exp(-r * t) - s * exp(-y * t) * pnorm(-d1)

    V
}

# Valuation with 1 stock price
put_option_pricer(s = 55, 60, .01, .02, 1, .05)
[1] 5.52021
# Valuation across multiple prices
put_option_pricer(s = 55:60, 60, .01, .02, 1, .05)
[1] 5.52021 4.58142 3.68485 2.85517 2.11883 1.49793

Let’s see what we can do with Rcpp. Besides explicitely stating the
types of the variables, not much has to change. We can even use the sugar function,
Rcpp::pnorm(), to keep the syntax as close to R as possible. Note how we are being
explicit about the symbols we import from the Rcpp namespace: the basic vector type, and
well the (vectorized) ‘Rcpp Sugar’ calls log() and pnorm() calls. Similarly, we use
sqrt() and exp() for the calls on an atomic double variables from the C++ Standard
Library. With these declarations the code itself is essentially identical to the R code
(apart of course from requiring both static types and trailing semicolons).

#include                                         
using Rcpp::NumericVector;
using Rcpp::log;
using Rcpp::pnorm;
using std::sqrt;
using std::log;

// [[Rcpp::export]]
NumericVector put_option_pricer_rcpp(NumericVector s, double k, double r, double y, double t, double sigma) {

    NumericVector d1 = (log(s / k) + (r - y + sigma * sigma / 2.0) * t) / (sigma * sqrt(t));
    NumericVector d2 = d1 - sigma * sqrt(t);
    
    NumericVector V = pnorm(-d2) * k * exp(-r * t) - s * exp(-y * t) * pnorm(-d1);
    return V;
}

We can call this from R as well:

# Valuation with 1 stock price
put_option_pricer_rcpp(s = 55, 60, .01, .02, 1, .05)
[1] 5.52021
# Valuation across multiple prices
put_option_pricer_rcpp(s = 55:60, 60, .01, .02, 1, .05)
[1] 5.52021 4.58142 3.68485 2.85517 2.11883 1.49793

Finally, let’s look at
RcppArmadillo. Armadillo has a
number of object types, including mat, colvec, and rowvec. Here, we just use
colvec to represent a column vector of prices. By default in Armadillo, * represents
matrix multiplication, and % is used for element wise multiplication. We need to make
this change to element wise multiplication in 1 place, but otherwise the changes are just
switching out the types and the sugar functions for Armadillo specific functions.

Note that the arma::normcdf() function is in the upcoming release of
RcppArmadillo, which is
0.8.400.0.0 at the time of writing and still in CRAN’s incoming. It also requires the
C++11 plugin.

#include 
// [[Rcpp::depends(RcppArmadillo)]]
// [[Rcpp::plugins(cpp11)]]

using arma::colvec;
using arma::log;
using arma::normcdf;
using std::sqrt;
using std::log;


// [[Rcpp::export]]
colvec put_option_pricer_arma(colvec s, double k, double r, double y, double t, double sigma) {
  
    colvec d1 = (log(s / k) + (r - y + sigma * sigma / 2.0) * t) / (sigma * sqrt(t));
    colvec d2 = d1 - sigma * sqrt(t);
    
    // Notice the use of % to represent element wise multiplication
    colvec V = normcdf(-d2) * k * exp(-r * t) - s * exp(-y * t) % normcdf(-d1); 

    return V;
}

Use from R:

# Valuation with 1 stock price
put_option_pricer_arma(s = 55, 60, .01, .02, 1, .05)
        [,1]
[1,] 5.52021
# Valuation across multiple prices
put_option_pricer_arma(s = 55:60, 60, .01, .02, 1, .05)
        [,1]
[1,] 5.52021
[2,] 4.58142
[3,] 3.68485
[4,] 2.85517
[5,] 2.11883
[6,] 1.49793

Finally, we can run a speed test to see which comes out on top.

s  matrix(seq(0, 100, by = .0001), ncol = 1)

rbenchmark::benchmark(R = put_option_pricer(s, 60, .01, .02, 1, .05),
                      Arma = put_option_pricer_arma(s, 60, .01, .02, 1, .05),
                      Rcpp = put_option_pricer_rcpp(s, 60, .01, .02, 1, .05), 
                      order = "relative", 
                      replications = 100)[,1:4]
  test replications elapsed relative
2 Arma          100   6.409    1.000
3 Rcpp          100   7.917    1.235
1    R          100   9.091    1.418

Interestingly, Armadillo comes out on top here on this (multi-core)
machine (as Armadillo uses OpenMP where available in newer versions). But the difference
is slender, and there is certainly variation in repeated runs. And the nicest thing about
all of this is that it shows off the “embarassment of riches” that we have in the R and
C++ ecosystem for multiple ways of solving the same problem.

To leave a comment for the author, please follow the link and comment on their blog: Rcpp Gallery.

R-bloggers.com offers daily e-mail updates about R news and tutorials on topics such as: Data science, Big Data, R jobs, visualization (ggplot2, Boxplots, maps, animation), programming (RStudio, Sweave, LaTeX, SQL, Eclipse, git, hadoop, Web Scraping) statistics (regression, PCA, time series, trading) and more…

Source:: R News

Leave a Reply

Your email address will not be published. Required fields are marked *

Time limit is exhausted. Please reload CAPTCHA.