Antonio Coppola edited C++ Functions.tex  over 9 years ago

Commit id: cade42599873fef26eeb64e1fc76bbb2a538061a

deletions | additions      

       

}  \end{lstlisting}  Aquick  microbenchmark comparison performed on OS X, Version 10.9.3, with a 2.9 GHx Intel Core i7 processor and 8 GB 1600 MHz DDR3 memory, reveals significant speed improvements: \begin{lstlisting}  microbenchmark(out.CPP <- lbfgs(objective(), gradient(), c(-1.2,1), invisible=1),  out.R <- lbfgs(objective.R, gradient.R, c(-1.2,1), invisible=1))  \end{lstlisting}  The results are the following: following, including also the optim routine as a benchmark:  \begin{lstlisting}  Unit: microseconds  expr min lq median uq max neval  lbfgs(objective(), gradient(), c(-1.2, 1), invisible = 1) 79.159 83.3620 88.1875 93.8975 189.992 100  lbfgs(objective.R, gradient.R, c(-1.2, 1), invisible = 1) 275.400 285.6185 303.1600 316.7580 1393.450 100  optim(c(-1.2, 1), objective.R, ..., method = "L-BFGS-B") 366.45 382.2585 415.9205 485.3085 20537.39 100  \end{lstlisting}  \subsection{Extra parameters in C++ implementations}  Much like in the R case, the passing of extra parameters with C++ implementations of the objective and gradient is achieved through the use of R environments. The following is an example replicating the Poisson regression using data from Golub et al. \cite{Golub_1999} with C++ function implementations. As before, we set up the objective and gradient as character strings. We include the extra environment argument and obtain data by evaluating it usign the \code{Rcpp::Environment} class:  \begin{lstlisting}  likelihood.include <- 'Rcpp::NumericVector lhood(SEXP xs, SEXP env){  arma::vec par = Rcpp::as(xs);  Rcpp::Environment e = Rcpp::as(env);  arma::mat X = Rcpp::as(e["X"]);  arma::vec y = Rcpp::as(e["y"]);  double prec = Rcpp::as(e["prec"]);  arma::mat Xbeta = X * par;  double sum1 = sum(y % Xbeta - log(1 + exp(Xbeta)));  arma::mat sum2 = sum(pow(par, 2 * prec));  arma::vec out = -(sum1 - 0.5 * sum2);  Rcpp::NumericVector ret = Rcpp::as(wrap(out));  return ret;  }  '  gradient.include <- 'Rcpp::NumericVector grad(SEXP xs, SEXP env){  arma::vec par = Rcpp::as(xs);  Rcpp::Environment e = Rcpp::as(env);  arma::mat X = Rcpp::as(e["X"]);  arma::vec y = Rcpp::as(e["y"]);  double prec = Rcpp::as(e["prec"]);  arma::vec p = 1 / (1 + exp(-(X * par)));  arma::vec grad = -((trans(X) * (y - p)) - par * prec);  Rcpp::NumericVector ret = Rcpp::as(wrap(grad));  return ret;  }'  \end{lstlisting}  \end{lstlisting}