Like many of us, I've written pipelines that implement best practices workflows for alignment of whole-genome sequencing data. These include a duplicate marking step, usually with Picard or Samtools aimed at identifying PCR and optical duplicates in raw sequencing data.
Illumina has a standard sequence identifier that are used for read names:
@[instrument_name]:[flowcell_name]:[tile]:[x]:[y]#[index]/[read_number]
These names, particularly the x and y coordinates, vital to detect optical duplicate reads.
The Sequence Read Archive (SRA) is an extremely valuable resource that has enabled genomics research to a staggering degree. The SRA also completely discards read names, a decision that is fairly controversial on the SRA github. Apparently you can recover original read names using the cloud data delivery service, but it costs money to download the data from s3/gcp, which adds up very quickly when dealing with sequencing data. Furthermore, this bit confuses me, because if the NCBI can deliver data with original read names, then clearly they are storing them, so why doesn't the standard SRA toolkit return them? But I digress.
I'm not going to argue with their decision here - that's been done plenty in other areas. There does appear to be little readily available data or discussion on the specific, actual impact of optical duplicates that are now impossible (or at least much more difficult) to identify once the data have been churned through the NCBI/SRA. I'm sure that such data exists, but either it's difficult to find or I'm very bad at finding it.
Are optical duplicates that remain in SRA data meaningfully harmful to downstream analyses? Can we do anything about them?
The -F option (or any equivalent) was removed in
fasterq-dump, andfastq-dumpis deprecated.Right, the intent of making this post is to discuss how much bioinformaticians should care about specifically optical duplicates. Technical details of sequencing runs, such as specific characteristics of duplicates, aren't always rigorously taught.
One thing I'm further curious about based on your reply, is because we distinguish PCR and optical duplicates, are the two meaningfully different in how they present in terms of their sequences? ie. if the reads are identical, then a script (such as
clumpify.sh) would be able to identify/remove both sets in one pass. Unfortunately in cases where one would want to _only_ remove optical duplicates, this wouldn't work, but in e.g. WGS this would be fine (and I'd expect standard dupmarking tools to remove the optical dups anyways, if this were true).Curios if there a reference for this?
fasterq-dumpseems to allow for parallel streams butfastq-dumpis still included in the sratoolkit. Fasterq-dump appears to have--use-nameand--seq-deflineoptions perhaps they would be relevant to recover the original headers (if present).Not in terms of sequence (unless the data has UMI's) but
clumpify.shcan distinguish between the two by using the cluster distance option.I am not sure how how big of a problem optical duplicates are. If sequencing facility is following loading recommendations/good practices they should not be a big problem. Do you have a specific reason/experience to be worried about them?
Apparently it's not deprecated yet, but will be soon and the SRA intends for it to be at some point "soon". In either case with some SRA data I've look at, the read names were still stripped even with -F - I suppose the data was originally submitted with stripped read names, but that just makes the original question all the more relevant because even if SRA tracked read names then a lack of coordinate information for reads is still a real world problem we encounter.
The source code of BBMap/clumpify indicates that clumpify computes cluster distance based on X and Y coordinates from an illumina-formatted read name. Therefore if the read names are not following Illumina spec, clumpify.sh cannot compute distances and distinguish optical from PCR duplicates.
I am not sure either if they're in practice a big problem, which brings us back to my original question, which in hindsight I could have phrased more generally:
Although if they present the same as PCR duplicates when coordinate information is absent, then one can probably assume that they would be caught and removed by dup marking software ie. end bam is the same from either case, although it might take a bit longer because the sequences need to be checked. But this won't work in cases where you wouldn't want to mark all dups (ie. amplicon sequencing)
With public data one has to work with what is there. If the information about position is absent then only thing that could be done is to look at sequence information. That can be done without alignment (with clumpify) or after alignment (with picard).
If removing duplicates is your aim then either of above would work. Nothing can be done about optical (or really cluster) dups in that case.