みっし～の研究生活

« 2009年02月 | メイン | 2010年01月 »

2009年08月11日

FASTA + FASTA.QUAL = FASTQ

以下，次世代シークエンサー関連の自分用のメモ

FASTAファイルとFASTA.QUALファイルからFASTQを生成するPerl (BioPerl)スクリプト

http://www.nabble.com/How-to-input-fasta-and-qual-into--Bio::Seq::Quality--td18634349.html を参照したけど，エラーが出るので，修正。BioPerl（というかPerl）についてはシロウトなので，手こずった。 "->"がRubyでいうところのメソッド呼び出しの"."なのね。

使い方は（スクリプトをfastaqual.plとすると）

perl fastaqual.pl hoge.fasta > fuga.fastq

でhoge.fastaとhoge.fasta.qualからfuga.fastqを作ります。 このスクリプトはhoge.fastaとhoge.fasta.qual内のフラグメントの順番が一致していることを前提としていることに注意。

#!/usr/bin/perl -w

use warnings;
use strict;
use Bio::SeqIO;
use Bio::Seq::Quality;

my ($seq_infile,$qual_infile)  =(scalar @ARGV == 1)
                         ?($ARGV[0]    ,"$ARGV[0].qual")
                         :@ARGV;

#Create input objects for both a seq (fasta) and qual file

my $in_seq_obj =
  Bio::SeqIO->new( -file   => $seq_infile,
		   -format => 'fasta',
		 );

my $in_qual_obj =
  Bio::SeqIO->new( -file   => $qual_infile,
		   -format => 'qual',
		 );

my $out_fastq_obj =
  Bio::SeqIO->new( -format => 'fastq'
		 );

while (1){
  ## create objects for both a seq and its associated qual
  my $seq_obj  = $in_seq_obj->next_seq || last;
  my $qual_obj = $in_qual_obj->next_seq;

  #use seq and qual object methods feed info for new BSQ object
  my $bsq_obj =
    Bio::Seq::Quality->new( -seq  => $seq_obj->seq(),
			    -qual => $qual_obj->qual(),
			    -id => $seq_obj->id(),
			  );

  $out_fastq_obj->write_fastq($bsq_obj);

}

投稿者 hmishima : 11:30