みっし～の研究生活

« CentOS 5.2のXVNCサーバーのGLX対応 | メイン | BioRuby版 FASTA + FASTA.QUAL = FASTQ »

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 : 2009年08月11日 11:30