2012/07/18 (Wed) R: 無理矢理pythonと連携しても、そこそこ速い

以前、RをCと連携させることで、特にループ処理の速度を改善させることができることを書いた（記事へのリンク）。しかし、やっぱりCは難しい。さらに、RからCへ渡すことのできる引数の型に制限があるなど、書き方の工夫もそれなりに手強い。

そこで、手っ取り早くループ処理の速度を改善する方法として、pythonと連携させてはどうかと考え、やってみた。pythonとの連携については特別なメソッドが用意されているというわけではない。以下の手順で、強引に連携させる。

1. 特定の計算を行うpythonのスクリプトを書き、結果はファイルに出力させる。
2. Rにおいて、system() 関数を使ってそのpythonスクリプトを走らせ、結果をファイルから読み取る。

system() 関数は、外部プログラムを呼び出すための関数で、引数としてコマンド（文字列）を受け取る。例えば、

system("ls")
とすれば、作業中のディレクトリのファイルリストを出力する（windowsなら system("dir") だろうか）。

したがって、
system("python file.py arg1 arg2")
とすれば、R上からpythonスクリプトを走らせることができる。もちろん好きなだけ引数を渡すことができる。ただし、戻り値を受け取ることはできないので、結果の受け渡しにはファイルを経由する。人に仕事頼むのだが事情があって結果を直接会って話すことができないのに似ている。結果を記したメモを壁に貼っておいてもらうのだ。


処理内容は以前と同じくN*N回の加減の計算。R, Cと連携, pythonと連携の3パターンの計算速度を比較する。
処理に用いたファイル一式はこちら。


下のグラフでは、横軸にNを、縦軸に処理にかかる時間をプロットしている（両対数軸）。やはりCとの連携が圧倒的に速い。pythonも健闘していて、Rにやらせるよりはかなり速い。

具体的には、N=10000の場合（つまり、1億回の計算）、R では145秒かかるところを、pythonでは15秒、Cではなんと0.1秒未満で終わらせている。























回帰分析の結果：

* C *
log(time) = -20.30163 + 1.93688 * log(N)  

* python *
log(time) = -14.66641 + 1.88273 * log(N)  

* R *
log(time) = -13.99015 + 2.06207 * log(N)  

 

Wald testを行うと、切片についてはすべての組み合わせについて、「等しい」という仮説は5％の有意水準で棄却された（Rとpythonの差はp値0.04くらいで、ぎりぎりだが）。
 

傾きについては、pythonとCの間には有意な差は見られなかった（p=0.156）。
一方、Rと比較したとき、python, Cはそれぞれ傾きが有意に小さい（それぞれp=0.002, 0.000）。



結論：Cほどではないがpythonとの連携もループ処理の速度改善に役立つ。

 

2012/02/01 (Wed) R: Cと連携するとここまで違う

RとC言語の連携の仕方を最近知った。友人の情報によると、サンタフェの人たちが言うには、この連携機能こそがMatlabにはない魅力なのだそうだ。

普通の回帰分析をやっている分には別にCなんて触れない方がいい。デバッグが面倒くさいだけだ。しかし、ちょっと面倒なことをやりたい時には、もしかしたらCと連携することで劇的な速度改善が見込めるかもしれない。


連携の良いところは、乗り換えではないということだ。原則としてRで作業して、計算が大変なところだけCの力を借りれば良い。例えば、データを読み込んで整形するところまでをRでやって、最尤法でパラメータを推計するために、最適化問題を解く所だけCのメソッドを利用する、という感じだろうか。やったことがないから想像の域をでないが、1日色々実験したところ、できそうな気がする。


今日は、Rをはじめとするインタープリタ言語が苦手であると言われているループ処理において、どれくらい差が出るのか試してみた。

作業内容はこちら↓↓

N <- 100
out <- 0

for (i in 1:N) {

    for (j in 1:N) {

        out <- out + i - j

    }

}

N * N の行列を、左上から右へ進み、各セルにおいてそのセルの行番号を足し、列番号を引く。右端に行き着いたら次の行へ・・・。対象なので答えは必ずゼロになるからプログラムとしてはまったく無意味なのだが、ポイントはN^2回の計算をするためNが増えると計算量が急激にあがるということだ。

結果は以下の通り。tC, tRはそれぞれCのメソッドを利用した場合, 通常どおりRでやった場合にかかった時間だ。差は歴然。正直ここまで違うとは思ってなかった。
          N    tC      tR

  1000 0.003   1.766

  2000 0.004   7.513

  3000 0.009  17.313

  4000 0.017  30.836

  5000 0.039  49.249

  6000 0.038  68.721

  7000 0.051  92.286

  8000 0.066 115.140

  9000 0.084 146.736

 10000 0.102 190.155

グラフにするとこんな感じ。

これじゃぱっと見よく分からないから、両対数にしてみる。すると、面白いことにほぼ平行な直線の上に乗っているようだ。






















最小二乗法で推計してみると、
log(tC) = -17.88 + 1.68 log(N)
log(tR) = -13.28 + 2.01 log(N)
という結果。tRの方は、2重ループの場合、理論的にlog(N)の係数は2になるんだろう。tCの方は、標準誤差が0.13なのでlog(N)の係数が2であるという仮説は、一般的な有意水準で棄却される。Cのメソッドを使った方がサイズに対する弾力性がほんの少し低い（のかな？）。むしろ、大きな差の源泉は切片の方にあるようだ。

2012/01/31 (Tue) R: コンソールでやるのもそこまで悪くない

Linuxの1つの弱点は、RのGUIがないことだと思っていたんだけど、コンソール上でやるのも意外に悪くないみたい。非公式のLinux用GUI（JGR, RKWard）はよく止まったり挙動がおかしくなるし、R Commanderも使いにくいことこの上ない。比べて、コンソールでやる分には何事も少なくともスムーズだ。唯一、ソースの一部を走らせるのが困難な点だけが残念なところ。そこはコードの書き方で工夫するか。


コンソールでRスクリプトを走らせるには、

source("file.R")

とする。でも、これだとあえてprint やcatを使わない限り何も表示されなくて味気ない。そういう場合には、

source("file.R", echo = T)

とする。

2012/01/30 (Mon) R: Commenting out multiple rows

R does NOT have.  Just do the following

if (F) {

......

}



【Rで複数行をコメントアウトするには】

2012/01/27 (Fri) R: R で計量経済やるときの大雑把な情報

こちら。

http://cran.r-project.org/web/views/Econometrics.html