Re^2: awk 2 perl oneliner

Thanks a lot, your script is great. Your solution is also a lot faster than the awk oneliner. I modified it a bit to fit my real data.

#!/usr/bin/perl

# create names lookup table from first file
while (<>) {
    #get rid of newline in name because there might be 
    #something like >name1\stext\n in the second file
    chomp;
    $names{$_} = 1;
    last if eof;
}
$/ = "\n>";  # set input record separator

# scan second file
print ">";            # print first ">" that would be missing
while (<>) {
      
    print if /^>?([\w\d]+)/ && $names{$1};
}
[download]

The only prob I have now is, when applied to my real data a ">" sign is missing at the beginning (which I got around with a print before the second while loop) and one is too much at the end.

Comment on Re^2: awk 2 perl oneliner Download Code

Replies are listed 'Best First'.
Re^3: awk 2 perl oneliner by almut (Canon) on Mar 06, 2010 at 20:16 UTC
The only prob I have now is, when applied to my real data a ">" sign is missing at the beginning ... and one is too much at the end. The reason for this is the way the records are being split. The first record (and only the first) will have a leading ">", and all records except for the last will have a trailing ">". Now, if you happen to not print out the first record (because it doesn't match) the initial ">" will be missing from the output file. Similarly, if you happen to not print the last record, there will be a stray ">" from the previous record. The proper way to handle those corner cases would be to always remove any leading or trailing ">"s, and always add a new one for output: `... # scan second file while (<>) { s/^>//; # remove leading '>' s/>$//; # remove trailing '>' print ">$_" if /^([\w\d]+)/ && $names{$1}; # add ^ }` [download] (This way, you also no longer need to allow a leading ">" in the regex (the `>?` part). )	[reply] [d/l] [select]
Re^4: awk 2 perl oneliner by space_agent (Acolyte) on Mar 07, 2010 at 12:57 UTC
Thank you very many, almut ;-). Now I have an efficient way to extract sequences from a fasta file (doing bioinformatics stuff). It is more than 10times faster than the awk oneliner from my first post. Here is almuts code slighly adjusted. It can be used to extract sequences from a multifasta file. `usage extract-seq.pl file1 file2 cat file1 AV12349 J234054 . . cat file2 >AV12349 AAACCGTACCTAGAACAGTA ACAGTACA >J434594 ACCGTTAGTACGATACA . .` [download] file1 contains ID to extract, file2 contains sequences in multifasta format. output will be written to the file extracted.seqs `#!/usr/bin/perl $output = "extracted.seqs"; open (OUT,">$output") or die; # create names lookup table from first file while (<>) { chomp; $names{$_} = 1; last if eof; } $/ = "\n>"; # set input record separator # scan second file while (<>) { s/^>//; #remove leading ">" s/>$//; #remove trailing ">" if (/^([\w\d]+)/ && $names{$1}) { print OUT ">$_";} }` [download]	[reply] [d/l] [select]