2
$\begingroup$

I tried to filter out FASTQ reads which are shorter than 259 bp with

bioawk -cfastx 'length() >= 259 {print "@" " " "\n""\n+\n"}' good/SZ005_NoIndex_L002_R1_009.good.fq.gz \
  | gzip >  good-filtered/SZ005_NoIndex_L002_R1_009.good-filtered.fq.gz`. 

Unfortunately, the output was wrong:

@ 

+

@ 

+

@ 

What did I miss?

Thank you in advance,

$\endgroup$
0
$\begingroup$

If you have paired-end reads, this solution keeps the two files in-sync (i.e. discard pairs where one of the two reads is shorter than 259). Also, it uses only Unix tools without requiring external programs:

paste <(zcat S1_L001_R1_001.fastq.gz | paste - - - -) <(zcat S1_L001_R2_001.fastq.gz | paste - - - -) \
| awk -v FS='\t' 'length($2) >= 259 && length($6) >= 259' \
| tee >(cut -f 1-4 | tr '\t' '\n' | gzip > r1.fq.gz) \
| cut -f 5-8 | tr '\t' '\n' | gzip > r2.fq.gz
$\endgroup$
2
$\begingroup$

You need to use variable $seq. For example:

bioawk -cfastx 'length($seq)>=259{print "@"$name"\n"$seq"\n+\n"$qual}' test.fq.gz

There are also more convenient tools like seqtk and seqkit. With those, you may

seqtk seq -L260 test.fq.gz
$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.