gcskewodysseyjapanese

Introduction

自分の解析の為に作ったサブルーチンであっても、G-language GAE では簡単にこれをG-language 標準解析関数に変換することができます。重複した解析ソフト開発はバイオインフォマティックスの発展を阻害するものであり、さまざまな解析手法を一つの場に集めるという意味でもG-language GAEの関数として解析手法を蓄積することには大きな意義があります。もちろん、あなたが作った解析手法を世界中の研究者が使うという可能性だけでも新鮮であることはいうまでもありません。 そこで、ここでは先程作成したGC skew のサブルーチンを標準関数化してみることにしましょう。

Step 0 標準関数の決まりを知る

多くの人が使うことになる標準関数には当然ながらいくつか決まりごとがあります。しかし、その多くは制約というよりは可能性を広げるものです。 詳しくは公式ドキュメント

標準関数 Manual、

SubOpt API Manual、

Messenger API Manual

を参照していただく必要がありますが、ここでは要約して説明します。

１、SubOpt APIを使いオプションを有効にする ２、オプションはデフォルト値を持たせ省略できるようにし、G のインスタンスの入力のみで使えるようにする ３、オプションは "output" と "filename" を含むのが望ましい ４、全ての出力はG::Messenger APIを使う ５、一つの関数は一つのサブルーチンから構成されるのが望ましい ６、グラフ作成には G::Tools::Graph::grapher を使用するのが望ましい ７、出来る限り簡潔に、みやすく、汎用的に、美しいコーディングをするのが望ましい

5 に関しては既に一つのサブルーチンになっており、6に関しては既に使用しており、7に関しては use strict; して努力するだけなので、ここでは説明しないこととします。

Step 1 SubOpt APIを使いオプションを有効にする

SubOptAPI とは、標準関数を使うときに -output=>;show" のように ハイフンから始まって矢印でつなぐ、GetOptに似たサブルーチン用オプションを可能にするAPIです。このAPIがG-language GAEの汎用性を高め、さらにシステム全体の統一性を実現しています。 SubOptの仕様は非常に単純明解です。opt_default()でデフォルト値を決定し、opt_get() でオプションを取得。さらにopt_val()でオプションの値を参照できるようになっています。

つまり、gcskew の関数ではウインドウがオプションですから、 &gcskew($gb, -window=>"1000"); というように使えたらいいわけですが、これは

  sub gcskew{ 
   my $gb = shift; 
   my $window = shift; 

を

  sub gcskew{ 
    my @args = opt_get(@_); 
    my $gb = shift @args; 
    my $window = opt_val("window"); 

のように変えてやることで可能になります。

まずopt_get()でオプションを取り出し、オプション以外の必要な引数を@argsにいれます。ここからGインスタンスをとりだし、opt_val("window") で別にオプションの値を取り出します。オプションの値は最初に取り出しローカル変数に保存しておくと便利です。

また、SubOpt APIにはGインスタンスに対して強力なメソッドが用意してあります。例えばこのGC skewの関数の場合、Gインスタンスで使用する部分は配列情報、つまり$gb->{SEQ} のみです。すると、ゲノム全体ではなく一部の配列のGC skewを計算したければその部分の配列をそのままいれたり、もしくは大きな配列の場合はそのリファレンスを代入することでメモリ効率を良く保ちたい場合があると思います。こんな時には

  sub gcskew{ 
    my @args = opt_get(@_); 
    my $gb = opt_as_gb(shift @args); 
    my $window = opt_val("window"); 

のように、opt_as_gb() 関数を使ってやると、それがGインスタンスでも、スカラー塩基配列でも、もしくはリファレンスでも自動的に$gb->{SEQ} を使えばその配列情報にアクセスできるように変換してくれます。

Step2 オプションはデフォルト値を持たせ省略できるようにし、G のインスタンスの入力のみで使えるようにする

さて、上のままでは -window オプションが入力されないと、Window が 0になってしまいエラーになってしまいます。だからといって、毎回window オプションを入力しないといけないのでは汎用的だとはいえません。そこで、opt_default() を使って初期値を入力してやります。

  sub gcskew{ 
    opt_default("window"=>10000); 
    my @args = opt_get(@_); 
    my $gb = opt_as_gb(shift @args); 
    my $window = opt_val("window"); 

こうしてやることで、windowオプションが入力されればその値が、入力されなければ 10000bp がウインドウサイズとして使用されます。

Step3 オプションは "output" と "filename" を含むのが望ましい

前回作った gcskew はファイルとグラフの二つの出力がありました。G-language GAE では容易に使えることに重点をおいていますので、標準関数関数には標準で全て output というオプションがついていて、"f" の値で'data'フォルダ内にファイル出力、"g"の値で'graph'フォルダ内にグラフ出力、指定なし、もしくは"show" の値でグラフを自動的に表示する、という機能を持つ必要があります。また、この時任意に出力するファイルの名前を指定できる必要があります。これにあわせて少しプログラム全体を変えてみましょう。outputの値次第で出力先を変える処理を加えてみます。 具体的には、opt_defaultの追加、opt_valをローカル変数に代入、そしてファイル名を$filenameに置き換え、出力先による設定を変えます。

  sub gcskew{ 
     opt_default("window"=>10000, "output"=>"show", "filename"=>"gcskew.png"); 
     my @args = opt_get(@_); 
     my $gb = opt_as_gb(shift @args); 
     my $window = opt_val("window"); 
     my $output = opt_val("output"); 
     my $filename = opt_val("filename");  
     my @gcskew = (); 
     my @location = ();  
     
     my $i = 0;   
     for ($i = 0; $i * $window < length($gb->{SEQ}); $i ++){      
           my $sequence = substr($gb->{SEQ}, $i * $window, $window);   
           my $c = $sequence =~ tr/c/c/;  
           my $g = $sequence =~ tr/g/g/;   
           my $skew = ($c-$g)/($c+$g);    
           push (@location, $i * $window);  
           push (@gcskew, $skew);      
     } 
     
     if ($output eq 'f'){ 
           mkdir ('data', 0777); 
           $filename =~ s/¥.png$/¥.csv$/; 
           my $j = 0; 
           open(OUT, '>data/' . $filename); 
           print OUT "location,GC skew¥n"; 
           for ($j = 0; $j <= $i; $j++){ 
                print OUT $location[$j], ",", $gcskew[$j], "¥n"; 
           } 
           close(OUT); 
     }elsif ($output eq 'g' || $output eq 'show'){ 
           mkdir ('graph', 0777); 
           G::Tools::Graph::grapher(¥@location, ¥@gcskew, -x=>'bp', 
                                           -y=>'GC skew', -title=>'GC skew',        
                                           -filename=>'graph/gcskew.png'); 
     
           system("gimv graph/gcskew.png"); 
     } 
         
     return @gcskew; 
  } 

Step 4. 全ての出力はG::Messenger APIを使う

だいぶここまででも標準関数らしくなってきましたが、もう一つ重要な決まりがあります。それは、全ての出力をG::Messenger APIを介してやらなければならない、ということです。G-language GAEは多彩なインタフェースを持ちます。プログラム言語としてコンパイラに読ませたり、コマンドラインで実行したり、グラフィカルユーザインタフェースを持ったり、さらにはウェブアプリケーションとしてさえ存在します。しかし、多彩なインタフェースでは全て同じ標準関数のプログラムが動いています。これを可能にするのがG::Messenger APIです。メソッドがMessengerに出力を渡すと、MessengerはG-language GAEが現在どのインタフェースを使って動いているかを判断し、的確な場所に出力を渡します。 Messenger で主に使われる関数は次の通りです。

msg_send() --- print STDOUT と同義 通常のアウトプットに。 msg_error() --- print STDERR と同義 エラーやシステムメッセージのアウトプットに。 msg_gimv() --- system("gimv "); と同義 グラフの表示に。 これでgcskew の標準関数化が完了です。

  sub gcskew{ 
     opt_default("window"=>10000, "output"=>"show", "filename"=>"gcskew.png"); 
     my @args = opt_get(@_); 
     my $gb = opt_as_gb(shift @args); 
     my $window = opt_val("window"); 
     my $output = opt_val("output"); 
     my $filename = opt_val("filename"); 
     
     my @gcskew = (); 
     my @location = ();  
      
     my $i = 0;   
     for ($i = 0; $i * $window < length($gb->{SEQ}); $i ++){      
           my $sequence = substr($gb->{SEQ}, $i * $window, $window);   
           my $c = $sequence =~ tr/c/c/;  
           my $g = $sequence =~ tr/g/g/;   
           my $skew = ($c-$g)/($c+$g);    
           push (@location, $i * $window);  
           push (@gcskew, $skew);      
     } 
     
     if ($output eq 'f'){ 
           mkdir ('data', 0777); 
           $filename =~ s/¥.png$/¥.csv$/; 
           my $j = 0; 
           open(OUT, '>data/' . $filename); 
           print OUT "location,GC skew¥n"; 
           for ($j = 0; $j <= $i; $j++){ 
                print OUT $location[$j], ",", $gcskew[$j], "¥n"; 
           } 
           close(OUT); 
     }elsif ($output eq 'g' || $output eq 'show'){ 
           mkdir ('graph', 0777); 
           G::Tools::Graph::grapher(¥@location, ¥@gcskew, -x=>'bp', 
                                              -y=>'GC skew', -title=>'GC skew',        
                                             -filename=>'gcskew.png'); 
        
           msg_gimv("graph/gcskew.png"); 
     } 
         
     return @gcskew; 
  } 

Step 5. To the next level ...

今回GC skew関数の標準関数化を説明しましたが、実際にG-language GAEに実装されているgcskew() 関数は AT skewもみれるようになっていたりと、汎用性やコードの簡潔さでより徹底した改良がされています。一度ソースコードをみて勉強してみてください。また、マニュアルにはより詳しい標準関数化の仕様について記載されていますので、こちらもじっくりと読んでみてください。 今後自分で作った解析メソッドを標準関数化し、G-language GAEに実装する機会が来ると思います。作った標準関数メソッドは、glang-devel@lists.sourceforge.jp に、サブジェクトが関数名、本文がソースコード、という形でメールすることで投稿できます。バイオインフォマティックスの未来の為にあなたが貢献する日は明日かもしれませんね。