Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Script extract_duplicates returned NaN duplication rate #1

Open
NgocTNguyen opened this issue Apr 25, 2017 · 2 comments
Open

Script extract_duplicates returned NaN duplication rate #1

NgocTNguyen opened this issue Apr 25, 2017 · 2 comments

Comments

@NgocTNguyen
Copy link

I ran ./extract_duplicates to estimate PCR duplication rate on our data:
./extract_duplicates --bam SRR1585519.dedup.merged.sorted.bam --VCF NA19099.het.vcf --mmq 20 > SRR1585519.dedup.hetreads

however it returned:
PCR duplicates marked 0 total-reads 0 frac nan discarded 62152388
#clusters
total reads (PE=1) 0 unique-reads 0 duplicates:0, duplication rate nan

Using samtool rmdup and IGV, we observed a lot of duplicate reads in our data, but the script here seemed to not work on our files. What's the problem here?

@bansal-lab
Copy link

is the data single-end reads ? If so, please use the "--singlereads 1" option.

@NgocTNguyen
Copy link
Author

NgocTNguyen commented Apr 27, 2017

I did add "--singlereads 1" in our command line, but the output was still the same.
I ran the script on another input file (PE reads) and did get the output file. So there is an issue of the script to take a SE file as an input.
However, when I ran the next script ./estimate_PCRduprate.py I got this output:

identified 0 het variants filtered 0 high-cov 1e-07 0.0 0.0 0.0
Traceback (most recent call last):
File "estimate_PCRduprate.py", line 407, in
[Frate,Ftotal,Funique,Funiqlist] = calculate_values(duplicate_counts,cluster_counts);
File "estimate_PCRduprate.py", line 103, in calculate_values
p = float(AB_counts[0])/(AB_counts[0] + AB_counts[1]);
ZeroDivisionError: float division by zero

How could it be "identified 0 het variants" ?

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

No branches or pull requests

2 participants