Which should not happen for a file for one sample (unless multiple sequencers are involved) :)
Hi, How can we do Read count from a batch of fastq.gz files.
Obviously if we unzip it then it is easy by grep.
But what is the best way to do without unzipping.
I tried this for a batch of files:
for i in `ls RIF*.fastq.gz | sort`; do zcat $i | echo $((`wc -l`/4)) ; done
Gives me result:
133408
133408
127328
127328
124862
124862
156815
156815
146333
146333
123459
123459
159771
159771
126039
126039
167112
167112
97320
97320
But When I unzip those same files and run grep with the identifier like
for i in `ls *.fastq | sort`; do grep "@M03691" $i |wc -l; done
I get true result as
132637
132637
126435
126435
123616
123616
156347
156347
145592
145592
123059
123059
158654
158654
125626
125626
166377
166377
97159
97159
Which is clearly a mismatch. Can anybody please suggest me what i am missing out in my first try. Thanks, Mitra
1 answer
Apparently not all of your read names start with M03691.
So far my money is on that as being the most likely explanation.
No ..thats not the fact ..and yes they all named as M03691 as much I can see.
Either they don't all start with that or you have cases where each entry isn't 4 lines. There are no other possibilities.
If this dataset was trimmed on the sequencer then there may be some entries with no sequence. Just a guess.
BTW, wc -l / 4 ends up truncating the division. It'd be interesting to not round and see if your files have line numbers that aren't multiples of 4.
Also when I am checking with one file I get exact match..but somehow above command not working..
grep "@M" RIFA3D26_S20_L001_R2_001.fastq |wc -l
145592
grep "@M0369" RIFA3D26_S20_L001_R2_001.fastq |wc -l
145592
This confirms that all the reads have same name.
Also with line count for single file matches.
cat RIFA3D26_S20_L001_R2_001.fastq | echo $((`wc -l`/4))
145592
Just still leaves me in confusion why line count don't work within zip files.
Your question isn't why line counting doesn't work (it does), but rather why you get different numbers. That's a rather important distinction.
Yes I described above two different ways of Read count from a batch of fastq.gz files ..out of which $((wc -l/4))
don't work (rather I must say not getting expected or correct result) when applied to zip files. So I assume this is something to do with the zip files ..and try to understand how the wc works in zip files.
Thanks Mitra
Log in to answer this question.
If you have
zcatthen you must havezgrep. Usezgrep -c "^@Sequencer_ID" file.fq.gzto get your counts.Great works thanks..don't know why zgrep didn't come to my mind.
But still I really wonder why wc -l`/4 didn't work. Definitely all these fastq files have 4 lines.