您好,登錄后才能下訂單哦!
這篇文章主要介紹“perl如何去除fasta或fastq文件中ID重復的序列”的相關知識,小編通過實際案例向大家展示操作過程,操作方法簡單快捷,實用性強,希望這篇“perl如何去除fasta或fastq文件中ID重復的序列”文章能幫助大家解決問題。
fastq和fasta 文件中序列ID都是唯一的,如果出現不唯一的情況,就需要給他去重復。這有一腳本可 實現此功能。
腳本幫助:
Usage Forced parameter: -fa fasta file <infile> -fq fastq file <infile> -fq1 fastq read1 file <infile> -fq2 fastq read2 file <infile> -f input file type,"fa"or"fq" <infile> must be given -od output dir <outfile> must be given -n output filename <outfile> must be given Other parameter: -h Help document
腳本既可以輸入fasta格式的,也可以輸入fastq格式序列文件。
-fa選項后跟輸入的fasta文件,-fq選項后跟輸入的fastq文件,如果fastq序列為雙端測序的,則-fq1后跟read1序列,-fq2后跟read2序列。-n后跟輸出文件前綴名,-od后跟輸出目錄。
腳本如下:
use Getopt::Long; use strict; use Bio::SeqIO; use Bio::Seq; #get opts my %opts; GetOptions(\%opts, "fa=s", "fq=s", "fq1=s", "fq2=s", "f=s", "od=s", "n=s", "h"); if(!defined($opts{f}) || !defined($opts{od}) || !defined($opts{n}) || defined($opts{h})) { print <<"Usage End."; Usage Forced parameter: -fa fasta file <infile>must be given -fq fastq file <infile>must be given -fq1fastq read1 file <infile>must be given -fq2fastq read2 file <infile>must be given -f input file type,"fa"or"fq" <infile>must be given -odoutput dir <outfile>must be given -noutput filename <outfile>must be given Other parameter: -h Help document Usage End. exit; } if($opts{f} eq "fa" && defined($opts{fa})){ my$read1 = $opts{fa}; open my $FQ1 ,"zcat $read1|" or die "$!"; my$fq1=Bio::SeqIO->new(-fh=>$FQ1,-format=>'fasta'); open my $GZ1 ,"| gzip >$opts{od}/${n}.fa.gz" or die $!; my$out1 = Bio::SeqIO->new(-fh => $GZ1 , -format => 'fasta'); my %id; while ( my $obj1=$fq1->next_seq() ) { my $id1=$obj1->id; if(exists $id{$id1}){ next; }else{ $id{$id1} = 1; } $out1->write_seq($obj1); } } if($opts{f} eq "fq"){ if(defined($opts{fq})){ my$read1 = $opts{fq}; open my $FQ1 ,"zcat $read1|" or die "$!"; my$fq1=Bio::SeqIO->new(-fh=>$FQ1,-format=>'fastq'); open my $GZ1 ,"| gzip >$opts{od}/${n}.fq.gz" or die $!; my$out1 = Bio::SeqIO->new(-fh => $GZ1 , -format => 'fastq'); my %id; while ( my $obj1=$fq1->next_seq()) { my $id1=$obj1->id; if(exists $id{$id1}){ next; }else{ $id{$id1} = 1; } $out1->write_seq($obj1); } }elsif(defined($opts{fq1}) && defined($opts{fq2})){ my$read1 = $opts{fq1}; my$read2 = $opts{fq2}; open my $FQ1 ,"zcat $read1|" or die "$!"; my$fq1=Bio::SeqIO->new(-fh=>$FQ1,-format=>'fastq'); open my $FQ2 ,"zcat $read2|" or die "$!"; my$fq2=Bio::SeqIO->new(-fh=>$FQ2,-format=>'fastq'); open my $GZ1 ,"| gzip >$opts{od}/${n}_R1.fq.gz" or die $!; my$out1 = Bio::SeqIO->new(-fh => $GZ1 , -format => 'fastq'); open my $GZ2 ,"| gzip >$opts{od}/${n}_R2.fq.gz" or die $!; my$out2 = Bio::SeqIO->new(-fh => $GZ2 , -format => 'fastq'); my %id; while ( my $obj1=$fq1->next_seq() and my $obj2=$fq2->next_seq() ) { my ($id1,$id2)=($obj1->id,$obj2->id); if(exists $id{$id1}){ next; }else{ $id{$id1} = 1; } $out1->write_seq($obj1); $out2->write_seq($obj2); } } }
關于“perl如何去除fasta或fastq文件中ID重復的序列”的內容就介紹到這里了,感謝大家的閱讀。如果想了解更多行業相關的知識,可以關注億速云行業資訊頻道,小編每天都會為大家更新不同的知識點。
免責聲明:本站發布的內容(圖片、視頻和文字)以原創、轉載和分享為主,文章觀點不代表本網站立場,如果涉及侵權請聯系站長郵箱:is@yisu.com進行舉報,并提供相關證據,一經查實,將立刻刪除涉嫌侵權內容。