I have PE125 sequence data and want to calculate the average read length of my reads.
samtools depth -q 0 -Q 0 my.bam | awk '{total+=\$3};END{print total}'The result show the average length of my reads is little smaller than 125, is this because that "samtools depth" can't calculate the inserted base?
If you just want to see average read length, this might help you
samtools view <in.bam> | awk '{ sum_length += length($10) } END { print sum_length/NR }'
If you are looking at insert length, Picard's CollectInsertSizeMetrics.jar might help you.
•
written
by
GouthamAtla
Thanks for your kind answer, I try your method to calculate the average length (result is 124.946, when I use samtools depth the result is 124.31), but it also not 125.What could cause this? And I also want to know if samtools depth can't calculate the inserted base?
Add: thanks geek_y, yes, we should consider clipping reads, it's the reason
You should take care of soft-clipping and hard-clipping.
Log in to answer this question.
