TSAに登録する際のmultifastaファイルの編集

NGSで読んだcDNAをデータベースに登録するとき、生配列はSRA（Sequence Read Archive）だけど、アセンブリ配列はTSA（Transcriptome Shutgun Assembly）という区分だそうな。

提出するfastaファイルは以下の形式に則っていないとダメだそうです。

 

>CLN01  <-- 1件目のエントリ名

ggacaggctgccgcaggagccaggccgggagcaggaagaggcttcgggggagccggagaa

ctgggccagatgcgcttcgtgggcgaagcctgaggaaaaagagagtgaggcaggagaatc

gcttgaaccccggaggcggaaccgcactccagcctgggcgacagagtgagactta

//      <-- 配列情報終了フラグ

>CLN02  <-- 2件目のエントリ名

ctcacacagatgcgcgcacaccagtggttgtaacagaagcctgaggtgcgctcgtggtca

gaagagggcatgcgcttcagtcgtgggcgaagcctgaggaaaaaatagtcattcatataa

atttgaacacacctgctgtggctgtaactctgagatgtgctaaataaaccctctt

//      <-- 配列情報終了フラグ

http://www.ddbj.nig.ac.jp/sub/mss/file_format-j.htmlより

perlってなにそれ状態で、大量の配列を含むfastaファイルの一括編集に手こずったので、メモ。

自分がやらないといけない作業は、スラッシュ//を追加することと、エントリ名を変更すること。次の2つのサイトを参考にコピペして、若干の変更を加えれば出来ました。なお、MacBook Pro・テキストエディタにmi使用。

マルチファスタの改行をとる(perl、awk) - script of bioinformatics

bioperl - Replace Fasta header files using Perl - Stack Overflow

使用したコードの一部は以下のもの。コピぺしたものをそのまま使っているので、コードの意味はイマイチわかってません。

#!/usr/bin/perl#!/usr/bin/perl

# マルチファスタファイル

open (IN, "<multi.fasta");

# 出力ファイル

open (OUT, ">output2.fasta");

@IN = <IN>;

$a = 1;

for ($i=0;$i<@IN;$i++) { 

  if(@IN[$i] =~ />/){

       if(@IN[$i] =~ /DN([0-9]*)/) {   #各配列のヘッダが"DNXXX"になっている

             $b = sprintf("%06d",$a);　#ヘッダの番号を6桁で表示"000001"のように

　　　  my $f = " / /\n\>\ABC$b\n";　#">ABC000001"という表記

             print OUT $f;

             $a =$a +1 ; 　　#番号を1ずつ足していく

       }else{

       print OUT @IN[$i];

       }

  }else{

print OUT @IN[$i];

}

}

close(IN);

close(OUT);