Fortranにやってもらう

Rをある程度習熟すれば、考えたことをなんでも実行できる。塩基配列の取扱いはBiostringsパッケージが強力だし、parallelパッケージの並列処理で多くの作業を高速化できる。

だけれども、やりたい計算があまりにもマニアックでapply系が使えないとき、とても困る。for()を多用せざるをえないような計算だ。

一応、for()も高速化はできる。先に適切なオブジェクトをつくっておき、その中のある座標だけを入れ替えるようにfor()ループを組み立てておけば、それなりに速く計算が終わる。宣言する感じ。そのつどオブジェクトを作り直すように組み立ててしまうと極端に遅い。それに、巨大なデータを用いて計算すると、速いといっても大したことない。

そこでFortranの出番となる。

RとFortranをつなぐのは.Fortran()という関数。Fortranのコードはプログラムではなくサブルーチンとして用意する。RからFortranに変数を渡せるし、FortranからRに変数を返してもらえる。ただし、サブルーチンの記述通りの順番で変数を渡す必要があるし、渡した変数しか返してもらえない。

サブルーチンの最後まで行った時点で保持していた変数のうち、Rから受け取ったものと同じもののみが、Rに帰ってくる。リストとして束ねられている。Fortranサイドで得た計算結果を受け取るためには、あらかじめRから初期値を入れた変数を与えて、それをFortranサイドで上書きしてしまえば良い。

RからFortranに変数が渡されると、Fortranが計算する。Fortranから変数が戻ってくるとき、Rのオブジェクトの中身だけを入れ替える形で戻ってくるらしい。だから、Fortranサイドでの宣言の仕方によってはおかしなことになる。

たとえば character(3) x みたいな感じで３個しか文字を格納できないようにしていると、Rから「x = "ABCDE"」というふうに５文字を受け取ったとしても Fortranで使うのは最初の３文字だけとなる（"ABC"）。この３文字を仮に別の文字に入れ替えてみるとする。たとえば、「x(3:3) = 'a'」とすると"ABa"になるのだが、このxをRに返してやると、上書きで返すのでRが受け取る変数は最初の３文字だけが改変された５文字、つまり"ABaDE"ということになる。