This is a test version of Biostars. For the public version, visit https://www.biostars.org.
extracting contigs

Hiya,

I have ordered my contigs by sequence length and wish to extract all of the contigs up to 310 (L50) or >000310F into a separate file.

Can this be done with grep?

Thanks!

fasta grep
• 1,581 views • link
updated by Ram Ram • written by hollyannj7 • hollyannj7

I'm not sure I understand whether you want the shorter or the longer than L50 contigs.

Anyhow, provided your contigs are ordered by contig length as you wrote, you can simply get the last line of the contig at L50 and head/tail the file, for example

head -n $LINE_L50 <contig_all.fna >contig_above_L50.fna
• 0 views • link
written by Carambakaracho Carambakaracho

Hi! Give a try to reformat.sh using maxlenght option (see BBMap reformat from BBMap toolkit).

• 0 views • link
written by iraun iraun

1 answer

The amazing seqkit can do this very easily:

seqkit seq seqs.fasta -m 310 > seqs.310.fasta

You can install seqkit using conda:

conda install -c bioconda seqkit

The docs are very good: https://bioinf.shenwei.me/seqkit/usage/#seq

• 0 views • link
written by Mark Mark

Log in to answer this question.

More posts like this