Hybrid-selection metrics from alfred qc

[evanbiederstedt]

Hi there

I was curious about how Alfred worked with hybrid-selection metrics. The motivation was we were discussing the differences between Alfred and Picard/GATK4's "CollectHSMetrics" for WES samples:

https://software.broadinstitute.org/gatk/documentation/tooldocs/4.0.5.1/picard_analysis_directed_CollectHsMetrics.php

java -jar picard.jar CollectHsMetrics \
      I=input_reads.bam \
      O=output_hs_metrics.txt \
      R=reference.fasta \
      BAIT_INTERVALS=bait.interval_list \
      TARGET_INTERVALS=target.interval_list

This requires both target and bait intervals.

It's not clear to what extent Alfred provides HS metrics, or whether the metrics output by e.g. CollectHSMetrics would differ that much from Alfred.

Somewhat related to this issue: #17

Any insight on your part is appreciated. Thank you for the help

[tobiasrausch]

I am not using picard but I would assume the bait intervals are the probes and the targets something like exons? Alfred only needs a single bed file of target regions.

[okbys]

I think there is a difference between Alfred and Picard when the mean depth is greater than 65535(0xFFFF).
The limit of "#AlignedBasesInBed" seems to be "#TotalBedBp * 65535".
Is there limit for "#AlignedBasesInBed"?

[tobiasrausch]

The basepair level coverage uses a 16bit integer by default. You can change this code fragment to uint32_t and recompile to allow coverages up to 2^32.

[okbys]

@tobiasrausch
Thank you for your solution.
However it's more complicated than I expected and it's not working...

At first, I changed this code from uint16_t to uint32_t, but alfred unexpectedly stopped without any output files.
Next, I changed this code from uint32_t to uint64_t, so alfred is never finished. (infinite loop?)
I checked the output file during the generation process, and the data in the ME section was as expected.

[tobiasrausch]

Let me check, I probably forgot somewhere in the code to refer back to this type definition.

[tobiasrausch]

You can specify now at compile time if you want to use 32bit integers for ultra-deep coverage. Just pull the last code and then use:

make all COVTYPE=uint32_t

Can you please try again if that works for you? Thanks.

[okbys]

It worked excellently!
The result is also as expected.
Thank you for your quick response!