Question: How to interpret the "Matrics" file of MarkDuplicates (Picard)
0
gravatar for ravipatel4
15 months ago by
ravipatel450
ravipatel450 wrote:

Hi,

I am using Picard to remove duplicates from BAM files. The manual says that it matters how the BAM files are sorted and you would get different results for different ways of sorting. I wanted to test how it affects with my data. So, I sorted my paired-end BAM file using either query name or coordinates using samtools (with -n option or without -n option, respectively). I get different output BAM files (different number of alignments) as mentioned in the manual. However, the duplication statistics in "Duplication Matrics" (output from "M" option) file is same for both. Shouldn't it be different for different ways the BAM files was sorted? Am I missing something? Please let me know. I have provided both commands below for your review.

Thank you in advance.

Best, Ravi

Command used for coordinate sorted BAM file:

time java -jar picard.jar MarkDuplicates I=coord.sorted.bam O=coordSorted_marked.bam M=coordSorted_marked.txt TAG_DUPLICATE_SET_MEMBERS=true TAGGING_POLICY=All REMOVE_DUPLICATES=true ASSUME_SORT_ORDER=coordinate READ_NAME_REGEX=null

Command used for query name sorted BAM file:

time java -jar picard.jar MarkDuplicates I= readName.sorted.bam O=querySorted_marked.bam M=querySorted_marked.txt TAG_DUPLICATE_SET_MEMBERS=true TAGGING_POLICY=All REMOVE_DUPLICATES=true ASSUME_SORT_ORDER=queryname READ_NAME_REGEX=null

Duplicate statistics for coordinate sorted BAM file:

Unknown Library 390015  32881469    0   6736693 292349  8862584 0   0.272361    49457941

Duplicate statistics for query name sorted BAM file:

Unknown Library 390015  32881469    0   6736693 292349  8862584 0   0.272361    49457941
markduplicates picard • 949 views
ADD COMMENTlink modified 15 months ago • written 15 months ago by ravipatel450

Picard markduplicates is marking reads as duplicate based on their sequence. Why should the result change whether you sort the file based on co-ordinates or name?

If you are looking to identify reads as optical duplicates (as opposed to PCR duplicates) based on their cluster co-ordinates, then you need to use clumpify.sh from BBMap suite: A: Introducing Clumpify: Create 30% Smaller, Faster Gzipped Fastq Files

ADD REPLYlink modified 15 months ago • written 15 months ago by genomax75k

Hi @genomax,

Because the Picard manual says so. My understanding of this based on their manual is as follows. In case of single-end reads, you can predict PCR duplicates by finding reads that map at a same location on the genome. In case of paired-end reads, you could use the paired-end information to more accurately predict PCR duplicates (by looking for read-pairs that have both ends mapped at a same location on the genome). But it might be computationally challenging to find read-pairs from a bam file unless the bam file is sorted by read names (reads in a pair have same name). My understanding is that Picard doesn't look for read pairs in the input bam file explicitly, unless they are one after the other. If the bam file is sorted by read name, the reads of a pairs will be one after the other, and Picard does the right thing of using paired-end information for calling PCR duplicates. If the bam file is sorted by chromosomal location, Picard may not use the paired-end information. Hence, it would matter how the bam file is sorted in case of paired-end data (which is the case in my data). I get different bam outputs (after removing PCR duplicates) depending on how I sorted the input bam file. However the matrix file is identical for different ways of sorting the input bam file. Therefore, I am trying to understand why is this the case.

Thank you for mentioning about clumpify.sh. I will keep that in mind if I wanted to find optical duplicates (Picard does that too by the way), but I am mostly interested in PCR duplicates at this point.

Thanks.

Best, Ravi

ADD REPLYlink written 15 months ago by ravipatel450

If an answer was helpful you should upvote it, if the answer resolved your question you should mark it as accepted.
Upvote|Bookmark|Accept

ADD REPLYlink written 15 months ago by Kevin Blighe52k
2
gravatar for h.mon
15 months ago by
h.mon29k
Brazil
h.mon29k wrote:

When marking duplicates in position sorted files, only mapped reads are analyzed. When marking duplicates in query name sorted files, unmapped reads will also be analyzed - unmapped reads may be marked as duplicates with query name sorted files. My interpretation of Picard MarkDuplicates (which may be incorrect) is that: 1) there is no difference in how mapped reads are analyzed for name or position sorted files, 2) for unmapped reads in name sorted files, probably only optical duplication is analyzed, but not PCR duplication.

So its not that there should be differences, but there may be differences, between the output of both methods. In your particular case, you disabled optical duplication detection (READ_NAME_REGEX=null), so I expect there will be no differences between both methods.

ADD COMMENTlink modified 15 months ago • written 15 months ago by h.mon29k

Hi @h.mon,

I see your point. But I think MarkDuplicates's behavior is different from your interpretation, because, there are different number of mapped reads in query name sorted input vs coordinate sorted input bam.

More importantly, the issue I am trying to resolve is the discrepancies between output bam files and output matrics files. My bam outputs are different (i.e. different number of output alignments) when I sort my input files in different ways. But the matrics file is same irrespective of the way I sort the input bam file. I am trying to resolve this discrepancy. I hope I could explain the issue clearly.

Thank you.

ADD REPLYlink modified 15 months ago • written 15 months ago by ravipatel450
2

You will see differences regarding exactly which reads/pairs are marked as duplicate depending on the sort order of the file, but the same number of reads/pairs will get marked as duplicates regardless (if that weren't the case the tool would be broken).

ADD REPLYlink written 15 months ago by Devon Ryan93k
1

How are the bam files different? How are you measuring or checking the differences?

but the matrix file is same

I believe you are talking about the metrics file, correct?

ADD REPLYlink written 15 months ago by h.mon29k
2

Hi h.mon,

Yes, I meant "matrics" file. Sorry about that.

I performed a test to see if the number of mapped reads are same in the output irrespective of the way the input is sorted. But that is not true. The MarkDuplicate's behavior seems different from your interpretation.

EDIT: I think I made a mistake in the command that used to count mapped reads in a bam file. It turns out that the number of mapped reads in output bam files is same regardless of how the file was sorted, meaning that your interpretation of how MarkDuplicates work looks correct. Sorry about the confusion.

ADD REPLYlink modified 15 months ago • written 15 months ago by ravipatel450
Please log in to add an answer.

Help
Access

Use of this site constitutes acceptance of our User Agreement and Privacy Policy.
Powered by Biostar version 2.3.0
Traffic: 1525 users visited in the last hour