This is a test version of Biostars. For the public version, visit https://www.biostars.org.
How To Grep Largest Contig From A Multi Fasta File

Hi everybody, can any one tell me please how to extract a largest contig from a multi-fasta file ?? using awk or grep ??

awk
• 13,990 views • link
updated by FadyNabil • FadyNabil • written by HG HG

2 answers

You can use Heng Li's bioawk and samtools:

Then the commands would be

# sort the sequences by length
$ cat w.fasta | bioawk -c fastx '{ print length($seq), $name }' | sort -k1,1rn | head -1
989    HR5V3UP02C00KT

# extract the sequence from the file 
$ samtools faidx w.fasta HR5V3UP02C00KT
>HR5V3UP02C00KT
TCGTACTCGTACGTAGAGGTTCGATCCTAGGGTCCTACGACGGAAGTAAAAACGGCCGGT
CCGGGCCCCGGTTCGACGTCGGACCGTAACCAACGAAAATTGGCCGGTAAAGGGGGTTCC
...
• 277 views • link
updated by Ram Ram • written by Istvan Albert Istvan Albert

Than you so much but could you please check the error

$ cat list4a.fasta | awk -c fastx '{ print length($seq), $name }'| sort -k1, 1rn | head -1
sort: invalid number after `,': invalid count at start of `'
• 0 views • link
updated by Ram Ram • written by HG HG

don't type in the whole thing at first, build it one step at a time, and pipe it through a pager, that way you will notice potential errors

• 54 views • link
written by Istvan Albert Istvan Albert

Dear Istvan, I am doing according to you suggestion:

cat list4a.fasta           #working fine 
awk -c fastx '{ print length($seq), $name }'

Here is the problem 

Usage: awk [POSIX or GNU style options] -f progfile [--] file ...
Usage: awk [POSIX or GNU style options] [--] 'program' file ...

could you please tell me waht is the error ??

• 0 views • link
updated by Ram Ram • written by HG HG

did you install bioawk? that's is the version of awk you will need to use.

• 387 views • link
written by Istvan Albert Istvan Albert

Hello Istvan, I have downloaded bioawk but when i am trying to install its showing 

hiren@FB11-10207:~/Desktop/spades/bioawk-master $ make 
make: `awk' is up to date.

But once i am checking bioawk like 

hiren@FB11-10207:~/Desktop/spades/bioawk-master $ bioawk
bioawk: command not found

Could you please let me know any solution ????

• 0 views • link
updated by Ram Ram • written by HG HG

can you explain this code please?

• 0 views • link
written by FadyNabil • FadyNabil 
cat seqs_oneline.fasta | perl -e 'while (<>) {$h=$_; $s=<>; $seqs{$h}=$s;} foreach $header (reverse sort {length($seqs{$a}) <=> length($seqs{$b})} keys %seqs) {print $header.$seqs{$header}}' | head -2

Source: sort fasta

• 0 views • link
updated by Ram Ram • written by Prakki Rama Prakki Rama

Log in to answer this question.

More posts like this