nameFile
seq_1014
seq_1039
seq_2848
seq_3213
seq_6847
seq_6980
seq_6997
seq_9319
seq_9561
seq_9850
# outputFile
>seq_1014
>seq_1039
>seq_2848
>seq_3213
>seq_6847
>seq_6980
>seq_6997
>seq_9319
>seq_9561
>seq_9850
>seq_10140
>seq_10141
>seq_10142
>seq_10143
>seq_10144
>seq_10145
>seq_10146
###Scriptuse
#!/usr/bin/env perl
my $list_file = $ARGV[0];
my $fasta_in = $ARGV[1];
my $fasta_out = $ARGV[2];
open(LIST_FILE, "<", $list_file) or die "could not open '$list_file' : $! \n";
open(FASTA_IN, "<", $fasta_in) or die "could not open '$fasta_in' : $! \n";
open(FASTA_OUT, ">", $fasta_out) or die "could not open $fasta_out : $! \n";
my @headers = ();
while(<LIST_FILE>) {
chomp;
next if ( /^\s*$/ );
push(@headers, $_);
}
my $pat = join '|', map quotemeta, @headers;
$/ = ">";
while(<FASTA_IN>) {
chomp;
if ( /$pat/ ) { print FASTA_OUT ">$_"; }
}
close(LIST_FILE);
close(FASTA_IN);
close(FASTA_OUT);
问题是不需要的输出
>seq_10140
>seq_10141
>seq_10142
>seq_10143
>seq_10144
>seq_10145
>seq_10146
我只想要确切的名称来匹配和 grep fasta 序列,但这个脚本很混乱,就像我只想要 seq_1014 并且它还给出了
>seq_10140
>seq_10141
>seq_10142
>seq_10143
>seq_10144
>seq_10145
>seq_10146
如何修复此脚本以获得所需的输出
答案1
那么两个问题:
- 不需要的匹配:这是因为您没有结束模式的每个子句以使其具有$最后表明您要匹配abced其次是没有什么。 seq_10140 做匹配seq_1014因为seq_1014 是在那里(某处)。你并没有说最后的另一个角色会导致失败seq_1014$。
- 当它应该匹配时失败:当我复制你的代码和输入文件时,这部分在我的系统上运行。 1039等等做匹配。我能想到的就是你$帕特由于您的模式文件存在一些问题,因此只有一个子句长,该问题正在通过将其剪切并粘贴到此论坛中或从该论坛中粘贴出来来修复。通常这意味着有两个行结束字符在每行末尾,就像 DOS 一样,整个文件在第一次读取时被读入。然而,它们还可能存在许多其他问题。将其进行调试或仅添加打印语句以查看其中的内容@标题和$帕特。通过“C”程序或“hexl-mode”下的 Emacs 获取文件的逐字节输出,看看其中是否有任何内容使您的读取感到困惑。