!!!課題2-2,2-3（課題2-1 発展編）
!!!0.1　課題2-2,2-3

*2-2：前章で作ったプログラムを改良し、担当遺伝子配列のアミノ酸配列を3通りすべて求める。

*2-3：相補的なDNA配列を求めるプログラムを書く。担当遺伝子の相補鎖のDNA配列を求め、そのアミノ酸配列を3通り求める。2-2の結果と合わせて6通りのReading Frameのうち、最も確からしい遺伝子領域はどれかを考察する。

!!!1　　課題2-2

みなさん先週の翻訳プログラムはうまくできたでしょうか。先週は配列ファイルの先頭から3つずつを順々にアミノ酸に翻訳していきました。

が、よく考えると、2番目の塩基から翻訳した場合、3番目の塩基から翻訳した場合では結果が異なってしまいます。

 C A T G C T G A C
 ~~~~~ ~~~~~ ~~~~~
 His    Ala   Asp
   ~~~~~ ~~~~~ ~~~~~
    Met   Leu   Thr
     ~~~~~~ ~~~~~ ~~~~~~
       Cys   Stop  

図のなかでは、第1のケースではこれは単なる配列です。しかし第2のケースではメチオニンが現れていますから、もしかしたらこれが遺伝子のはじまりかもしれません。さらに第3のケースではストップコドンが現れています。これが遺伝子の途中だとしたら、ここでコード領域が終わることになります。みなさんが切り取ったファイルはマイコプラズマ菌のゲノムを何の根拠もなくぶった切っただけですから、その先頭から3つずつ読んでいってきちんと遺伝子に探し当たる証拠はありません。

そこで、この3つのケースすべてを試してみる必要があります。

方法は、2つあります。

*先頭の塩基をいくつか削って残りの2通りファイルをつくり、前回つくったプログラムにかける

*3通りのケースすべてをプログラムで処理する

どちらの方法でもよいのですが、どうせですからここではプログラムを書くことにします。1回書いてしまえばあとは実行するだけですから楽です。

!!1.1　前回のプログラムをもとに翻訳

前回つくったなかで、実際に3文字ずつ翻訳するサブルーチンを作ったのですから、これを改良して塩基配列の1文字目からと、2文字目からと、3文字目からに対して実行すればいいわけです。

つまり、for文の開始がいままでは $i = 0 になっていたと思いますが、この初期値を$j などにしてやり、もう一つ外側にfor文を作り、$jが0から2までを繰り返してやれば3つの読み枠に対応できます。この時、返すアミノ酸配列も3つになるのですから、せっかくなのでこの$jを使って@aminoという配列に記憶させてみましょう。こうしておくと、最後のreturnで配列を返せばいいので楽です。

!!!2　　課題2-3

課題2-2でreading frameの3つの可能性をすべて表示することができるようになりました。しかし、実は可能性はまだあります。それは、DNAの二重螺旋の反対側にある場合です。

     5' __G T A C G A C T G __3'
     3'   C A T G C T G A C   5'

DNA鎖は分子内のリン酸の結合位置によって5'末端と3'末端の方向性があります。遺伝子は、その発現機構の性質で、5'側から3'側へと向かう方向にしかありません。

DNAの二重螺旋は2つの逆の方向を持った鎖によって構成され、互いに相補的なAとT、GとCによって互いに結合しています。

というわけで、いまある片方の配列からもう片方の配列を得るには、

1.いまある配列を逆向きにする

2.それから塩基をすべて相補的な塩基に置き換える

という操作をしなければなりません。
こうして得られた配列の3通りを加えて合計6通りの読み方を考慮に入れれば、完璧です。

!!2.1　相補鎖を求める関数

まず配列をひっくり返さなければなりませんが、これにはせっかくなので新しい関数を作ってしまいましょう。$seqを引数で渡すと、その相補鎖を返してくれる関数ができれば便利です。このために必要な手順は、配列の順序を入れ替えることと、AならばT、GならばCといった具合に相補塩基を当てはめることです。

まず、関数の枠組み自体は前回同様に作ります。

 sub complemental(){
        my $nuc = shift;
        my $complement = '';
 
        ************** #裏返す
        ************** #相補塩基を当てはめる
 
        return $complement;
 }

こんな時、文字列の順序を入れ替えてくれる関数があったら便利ですよね。しかし、そんなうまい話が、、、あるんです。Perlには標準でreverseという関数があり、$a = reverse("hoge"); としてやると$aに"egoh"が入ります。Larry Wallに感謝ですね。

さて、さらに相補塩基を当てはめる関数なんかがあればうれしいのですが、それはまだG-languageやBioPerlなどが必要で、標準では存在しません。しかし、Perlはもちろんこれを簡単にする仕組みを用意してくれています。そう、勘のいい人は気が付いたでしょう。tr///です。

tr///は一文字置換の関数なのですが、これは非常に便利で、特定の文字がきたら特定の文字に、という複数の条件を一文で書き表すことができます。つまり、

 $nuc =~ tr [atgc]
            [tacg];

と書くだけで、aをtに、gをcに置換してくれるのです。Perlはこのように非常に気の利いた言語ですので、とことん楽をしようと考えれば考えるほど、簡単に、綺麗なプログラムが書ける可能性が高くなると言っても過言ではないでしょう。

それでは、プログラムを書いてみましょう。できれば前回学んだように、できるだけ短くまとめてみましょう。

!!!3　　遺伝子（コード領域）を見つける

これでみなさんの担当の塩基配列について、全ての可能性でのアミノ酸翻訳結果が出揃ったことになります。では、生体内で実際にアミノ酸に翻訳される部分はどうやって見つければいいでしょうか。

厳密には、これはウエットな実験によって決定されるべきものですが、塩基配列の字面からだけでもある程度有効な予測をすることが可能です。

そのルールは以下のようなものです。

1.コード領域はatg（M,メチオニン）からはじまる。

2.バクテリアの場合はまれにgtgからはじまる。これは本来バリンのコードだが、スタートコドンとなる場合に限ってメチオニンに翻訳される。

前回、これをVと表示するようにプログラムしたはずです

3.コード領域はtaa,tga,tagいずれかのストップコドンが終端となる。

4.複数の可能性がある場合、1番長いものである可能性が高い。

ほとんどの場合、このルールに沿って探していけば、これだ！というものが見つかるはずです。これだけでは判断が難しい場合は、ケースバイケースでしょう。

これをプログラムで自動化してみてもいいですね。興味と余力のある人はやってみるとよいでしょう。

!!!4　　おまけ・より高度なプログラミングのために

遺伝子情報処理の扱う範囲はあくまでBiologyです。コンピュータをこの分野に使うことは決して単なる手段ではなく、コンピュータを使うことこそが本質的な事柄となり得るのですが、しかし実際の研究ではプログラムを書くこと自体は副次的な作業です。

とすれば、だからこそ、不必要に研究期間を延長することなく、また相互にプログラム資源を共有するために、読み易く、1回限りの使用で終わってしまわないプログラムを書くことが重要となります。

なるべく汎用性を意識し、繰り返す処理は関数化していくように心がけると良いでしょう。そして、ある程度Perlに慣れてきたら、「Effective Perl」を読むことを薦めます。きっと数倍効率よくPerlをかけるようになるはずです。

それでは、これからの研究を頑張ってください。