Rcppを触り始めた

Rcppを使うと従来の方法と比べてよりシームレスにRの拡張を作れるようになる.

RをC++で拡張する動機

主に以下の二点が挙げられる.

  1. 計算速度を向上させる
  2. C++から使えるライブラリをRで利用する

1については以下の記事が参考になる.

インストール方法

Rcppパッケージのインストール方法は以下の通り.

> install.packages("Rcpp")

インストール後,以下のコマンドでRcppの説明を見ることができる.

> vignette("Rcpp-introduction")

従来のC++によるRの拡張

ある決まった数値ベクトルを返す関数をC++で書き,Rで呼び出してみる.

numeric_vector.cc

従来のC++によるRの拡張の書き方.

#include <Rinternals.h>

extern "C" SEXP numeric_vector() {
  SEXP v; 
  PROTECT(v = allocVector(REALSXP, 2));
  REAL(v)[0] = 123.45;
  REAL(v)[1] = 67.89;
  UNPROTECT(1);
  return v; 
}

Rのガベージコレクションにオブジェクトを勝手に回収されないようにするために,自分で作ったオブジェクトをPROTECT/UNPROTECTする必要がある.この処理が非常に煩わしい.特にUNPROTECTについては,「ポインタを渡してオブジェクトを保護解除」するのではなく,「保護されている直近n個のオブジェクトを保護解除」するので,プログラマー自身がオブジェクトをPROTECTした順番をきちんと把握しておかなければならない.また,ベクトルに実数を代入する処理も決してきれいに書けるとは言えない.このような事実がC++でRの拡張を書く敷居を上げていると思う.

コンパイル

コンパイルはRコマンドを使って行う.

$ R CMD SHLIB numeric_vector.cc 

numeric_vector.soが生成される.

実行

Rのコード中でnumeric_vector.soをloadし,.Callで関数名を指定して呼び出す.

> dyn.load("numeric_vector.so")
> numeric_vector <- function() { .Call("numeric_vector") }
> numeric_vector()
[1] 123.45  67.89

RcppによるRの拡張

Rcppを使うと先程のソースは以下のように書き換えることができる.

numeric_vector_rcpp.cc
#include <Rcpp.h>

RcppExport SEXP numeric_vector_rcpp() {
  Rcpp::NumericVector v(2);
  v[0] = 123.45;
  v[1] = 67.89;
  return Rcpp::wrap(v);
}

RcppのコードではPROTECT/UNPROTECTは必要ない.これだけでかなりプログラムが書きやすくなるはずだ.また,Rcpp::NumericVectorはC++のvectorと同じように値を代入することができ,REAL(v)[0]のような煩わしい書き方も必要ない.

さらに,C++0xの初期化リストをサポートするコンパイラ(例:GNU g++ >= 4.4)では以下のように書くこともできる.

#include <Rcpp.h>

RcppExport SEXP numeric_vector_rcpp() {
  Rcpp::NumericVector v = {123.45, 67.89};
  return Rcpp::wrap(v);
}

素晴らしい.

あと,ここには示していないがRとC++の文字列(string)もかなりシームレスに扱うことができる.

Makevars

R CMD SHLIBでRcpp用のファイルをインクルード・リンクするため,あらかじめ以下のようなMakevarsファイルを用意しておく必要がある.

PKG_CXXFLAGS=$(shell Rscript -e "Rcpp:::CxxFlags()")
PKG_LIBS=$(shell Rscript -e "Rcpp:::LdFlags()")
コンパイル

先程と同様.

$ R CMD SHLIB numeric_vector_rcpp.cc
実行

これも先程と同様.

> dyn.load("numeric_vector_rcpp.so")
> numeric_vector_rcpp <- function() { .Call("numeric_vector_rcpp") }
> numeric_vector_rcpp()
[1] 123.45  67.89

以上

これだけの投資で大幅な速度向上が見込めるのであればやる価値は十分にあるだろう.

Tsukuba.R#8の発表ネタに使おうかな.