« 08.17 民主党の表は見にくい | ココ | 08.28 SVGってどうなのよ »
2009年8月19日
EMBOSSを試す その8
$ perl -lane 'BEGIN{$starts[0]=0;$offset=0;$margin=0};END{foreach(@res){$_->[2]-=$margin, print @{$_}}};exit if(/^#/);$margin=$-[1],$offset+=$+[1]-$-[1],push @starts,$offset if(/^ +(-[-:|]*)/); push @res,[$starts[$#starts/2]," ",$-[0]," ",$#starts%2," ",$&] while($_ =~ m/(?<= |\|)[A-Z][a-z]\w+/g);print' gq149630.remap
(↓パイプ直結ver.)
$ remap genbank:GQ149630 stdout -enzymes all -sitelen 4 | perl -lane 'BEGIN{$starts[0]=0;$offset=0;$margin=0};END{foreach(@res){$_->[2]-=$margin, print @{$_}}};exit if(/^#/);$margin=$-[1],$offset+=$+[1]-$-[1],push @starts,$offset if(/^ +(-[-:|]*)/); push @res,[$starts[$#starts/2]," ",$-[0]," ",$#starts%2," ",$&] while($_ =~ m/(?<= |\|)[A-Z][a-z]\w+/g);print'
ワンライナー考えてみた(長いけど)。その7で使ってみた制限酵素切断地図を出力するremap、きれいに整形して出力してくれるのはうれしいけど、そんなんよりも表として出して整理したいんだ! という人もいると思うんよね。gq149630.remapは、その7で出力したデフォルトの設定状態での結果(パイプ直結版と同じ)。もっといい解あったらぜひご連絡を。
上のプログラムを通すと、下のように「シーケンスのオフセット」「切断位置の直前(Fw)または切断位置の直前(Rv)」「Fw(=0)かRv(=1)か」「制限酵素名」の4項目が出てくる。同じシーケンスのオフセットの中に、同名の制限酵素が最大2回出てくるので、1回目が上の段(3項目目が0)、2回目が下の段(3項目目が1)のものということで。次回はもうちょっと「らしく」しようかね。
0 11 0 SpeI
0 12 0 MaeI
0 26 0 FaiI
0 36 0 CviRI
0 15 0 TspDTI
0 28 0 FaiI
0 40 0 AciI
0 48 0 CviRI
0 59 0 FaiI
0 15 1 SpeI
0 28 1 FaiI
0 36 1 CviRI
0 42 1 AciI
0 48 1 CviRI
0 59 1 FaiI
0 14 1 MaeI
0 26 1 FaiI
0 13 1 TspDTI
60 8 0 SetI
60 11 0 FatI
60 33 0 AlwNI
60 44 0 TatI
60 12 0 CviAII
60 34 0 SfeI
... (省略) ...
1620 37 1 Hpy178III
1620 20 1 PasI
1620 20 1 SecI
1620 19 1 BseBI
1620 19 1 ScrFI
1680 1 0 Tsp4CI
1680 6 0 TscAI
1680 11 0 FaiI
1680 17 0 MseI
1680 2 1 Eco31I
1680 11 1 FaiI
1680 19 1 MseI
1680 2 1 BsmAI
1680 1 1 SfeI
1680 0 1 Tsp4CI
By ただ at 21:06 カテゴリー ; Perl , 生命科学 , プログラミングとか
« 08.17 民主党の表は見にくい |
08月の記事
| 08.28 SVGってどうなのよ »
トラックバック
このエントリーのトラックバックURL:
http://pinmarch.sakura.ne.jp/mt/mt-tb.cgi/1463
