This is a test version of Biostars. For the public version, visit https://www.biostars.org.
extracting reads from fastq file based on read_id

I have a file with a set of read ids called tmp. I want to extract reads with these read ids from another fastq file called reads.fastq. When I use grep -f tmp reads.fastq I get only the first line (header) of the read. How can I get the output in full fastq format.

next-gen
• 8,726 views • link
updated by Charles Plessy Charles Plessy • written by sumithrasank75 sumithrasank75

3 answers

grep -A 3 seq_ID your file should get you the full fastq record. Iterate over the ID's you have.

filterbyname.sh from BBMap would also work with your list of ID's in a file.

• 1 views • link
written by GenoMax GenoMax

Today I learnt about -A. Thanks.

• 1 views • link
written by poisonAlien poisonAlien

If your tmp and reads.fastq have their reads in the same order, you can try our syncpairs tools.

• 1 views • link
written by Charles Plessy Charles Plessy

For bash, I would use a nice one-liner:

for x in $(cat tmp); do reads.fastq | grep $x -A 3; done

Of course you can redirect to wherever you want:

for x in $(cat tmp); do reads.fastq | grep $x -A 3; done > tmp2

I hope that this helps!

Brent Wilson, PhD | Project Scientist | Cofactor Genomics

4044 Clayton Ave. | St. Louis, MO 63110 | tel. 314.531.4647

Catch the latest from Cofactor on our blog.

• 1 views • link
written by brent_wilson brent_wilson

For m desired reads out of a fastq file with n total reads, this takes O(mn) time, as opposed to the -f flag given to grep (or just using BBMap), which is likely O(n).

• 1 views • link
written by Devon Ryan Devon Ryan

Log in to answer this question.

More posts like this