First, I define my NSE function in C++. This is done conveniently with the Rcpp package.
#include <Rcpp.h>
using namespace Rcpp;
//[[Rcpp::export]]
double NSE(NumericVector yhat, NumericVector y) {
// Rcpp comes with a lot of syntactic sugars so you can write C++ code almost like
// R code. Unfortunately it does not have the sugar '^' element-wise exponent yet, so
// we use x * x instead of x^2.
double ybar = mean(y);
double rss = sum((y - yhat) * (y - yhat));
double tss = sum((y - ybar) * (y - ybar));
return 1 - rss/tss;
}
Now let’s generate two random vectors to calculate the NSE.
x <- rnorm(100)
y <- rnorm(100)
NSE(x,y)
## [1] -1.020726
I’ll benchmark my NSE function with the one written in the hydroGOF package. First, let’s see whether the two functions produce the same answer.
hydroGOF::NSE(x, y)
## [1] -1.020726
Yes! Now, more importantly, let’s compare the speed.
microbenchmark::microbenchmark(hydroGOF::NSE(x, y), NSE(x, y), times = 1000)
## Unit: microseconds
## expr min lq mean median uq max neval
## hydroGOF::NSE(x, y) 19.7 20.8 24.1455 22.2 23.9 93.7 1000
## NSE(x, y) 1.3 1.6 3.2297 1.9 2.8 846.9 1000
The C++ implementation is about 15 times faster!
To learn more about using C++ in your code, please refer to the following sources: