MSK-ACCESS Quality Control (QC) Version 1 (v1) is intended to ensure all samples are processed in an accurate and consistent manner using the https://github.com/mskcc/ACCESS-Pipeline workflow.

The following pages give descriptions on how to interpret Quality Control metrics presented in the PDF report, page by page.

Here is an example Quality Control report for the V1 ACCESS pipeline:

Theoretical Method

Unlike other coverage metrics from this report which report coverage for bait regions, this graph shows the coverage of actual genomic target regions of the ACCESS A panel

Technical Methods

• Tool Used:

  • Marianas

  • Waltz CountReads

  • aggregate_bam_metrics.sh

  • tables_module.py

  • plots_module.r

• Input

  • Duplex Bams

  • pool A bed file

• Output

  • waltz_duplex_a_exon_level_files (directory of Pool A Exon Targets QC results)

  • waltz-coverage.txt

Interpretations

Coverage in this graph should be slightly higher than for the probe-level coverage results, as the calculation is limited to a smaller window of the histogram of coverage values. This metric is relevant for analysts who are more interested in coverage for a particular gene rather than coverage of the baits used to target that gene.

Theoretical Method

Total number of reads sequenced. Obtained from iterating through SAMRecord instances in Standard Bam file.

Technical Methods

• Tool Used

  • Waltz.jar CountReads

  • tables_module.py

Interpretations

ACCESS is designed to target ~50-80M reads per ctDNA sample, and ~5-10M reads for buffy coat samples. This number should be independent of the amount of library input DNA as well as coverage values, as PCR will bring low input values up to a consistent amount for sequencing.

Theoretical Method

This metric is obtained by iterating through the bam file, and looking at the sam flag which indicates whether each read has an adequate mapping to the HG19 reference.

Technical Methods

Waltz uses a method from the of the HTSJDK library:

• Tool Used

  • waltz.jar CountReads

  • Aggregate_bam_metrics.sh

Interpretations

Mapping fraction to the human genome should be above 97%, in most cases if it is below that, there is a chance that there is contamination from another species.

Theoretical Method

Divide the number of reads mapping to ACCESS genome bait regions by the sample’s total read count.

Technical Methods

• Tool Used:

  • waltz.jar CountReads

  • aggregate_bam_metrics.sh

Interpretations Ideally should be between 60%-80% and should not drop below 50% for ctDNA samples (A + B targets combined). Because A and B targets are mixed in 50:1 ratio, there should be a larger on-target rate for the A targets. If the rate drops below 50%, with adequate read counts, this could be indicative of a bad capture. For buffy coats, this is about 35%-45%.

Theoretical Method

Bin GC content of each region in the bam file into 5% intervals, and plot mean coverage across all regions that fall into each bin.

Technical Methods

• Tool Used:

  • Waltz CountReads

  • aggregate_bam_metrics.sh

  • tables_module.py

  • plots_module.r

• Input

  • Standard bam

  • Collapsed unfiltered bam

• Output

  • sample_id-intervals.txt

Interpretations Extreme base compositions, i.e., GC-poor or GC-rich sequences, lead to an uneven coverage or even no coverage of reads across the genome. This can affect downstream small variant and copy number calling. Both of which rely on consistent sequencing depth across all regions. Ideally this plot should be as flat as possible. The above example depicts a slight decrease in coverage at really high GC-rich regions, but is a good result for ACCESS.

Technical Method

• Tool Used: plots_module.r

Theoretical Method

Insert size is calculated from the start and stop positions of the reads after mapping to the reference genome.

Technical Methods

• Tool Used:

  • Waltz CountReads

  • aggregate_bam_metics.sh

  • tables_module.py

  • plots_module.r

• Input

  • Collapsed all unique bam

  • ACCESS pool A bed file

• Output

  • sample_id.bam.fragment-sizes

  • fragment_sizes.txt (aggregated across samples from a single bam type / pool combination)

Interpretations Cell free DNA has distinctive features due to the natural processes behind its fragmentation. One such feature is the set of 10-11 bp fluctuations that indicate the preferential splicing of fragments due to the number of bases per turn of the DNA helix, which causes a unique pattern of binding to the surface of histones.

The more pronounced peak at 166 bp indicate complete wrapping of the DNA around the histones’ circumference, and similarly the second more pronounced peak indicates two complete wraps.

Buffy coat samples are mechanically sheared and thus do not exhibit these distinctive features, hence the different shape for their distribution.

Theoretical Method

Calculate average coverage of each of four possible bam types for each sample:

• Standard or “Uncollapsed”

• Collapsed unfiltered: after merging all reads from same UMI family

• Collapsed simplex: three or more reads found on one strand

• Collapsed duplex: one or more reads found on both strands (top and bottom)

Coverage is first averaged across each position in a single bait region. Then, the average across each bait region in the sample represents the sample’s final coverage value.

Technical Methods

• Tool Used:

  • Marianas

  • Waltz CountReads

Interpretations

Expected range of coverage values:

• • Raw coverage A panel:

    • ctDNA: ~ 15000x-20000x

Theoretical Method

Similarly to the Pool A metrics, the UMI family type composition is here presented for Pool B targets. Buffy coat samples should have comparable UMI family composition for both Pools A and B.

Technical Methods

• Tool Used:

  • Marianas

  • make_umi_qc_tables.sh

  • plots_module.r

• Input

  • Marianas collapsed fastqs

• Output

  • family-types-B.txt

Interpretations

Duplex families are valuable for their low noise rate after collapsing, thus we'd like to see as high of a duplex "saturation" as possible. Because Pool B probes are mixed at a lower ratio in the capture process for cfDNA samples, they will have less duplex saturation. If this value is lower, we may not have captured enough of the original molecules to find both strands after PCR replication.

Theoretical Method

The sequencer uses the difference in intensity of the fluorescence of the bases to give an estimate of the quality of the base that has been read. The BaseQualityScoreRecalibration (BQSR) tool from GATK recalculates these values based on the empirical error rate of the reads themselves, which is a more accurate estimate of the original quality of the read.

Technical Methods

• Tool Used:

  • GATK BaseQualityScoreRecalibration

  • Picard MeanQualityByCycle

• Input

  • Standard, Uncollapsed Bams

• Output

  • sample_id.bam.quality_by_cycle_metrics

  • sample_id.bam.quality_by_cycle.pdf

Interpretations

It is normal to see a downwards trend in pre and post-recalibration base quality towards the ends of the reads. Average post-recalibration quality scores should be above 20. Spikes in quality may be indicative of a sequencer artifact.

Theoretical Method

Marianas performs read grouping based on the 6-base UMI sequence (three from each side of the DNA fragment), as well as the fragment start position (and stop position?). If multiple read pairs have the same information for these two metrics, they will be grouped into the same UMI "family".

UMI family types are defined by the following categories:

• Duplex: both top and bottom strand were found for this fragment

• Simplex: only one of (top|bottom) strand was sequenced, and >=3 copies for that strand were found

• Sub-Simplex: exactly 2 copies of a single strand were found

• Singletons: exactly 1 copy of a single strand was found

Technical Methods

• Tool Used:

  • Marianas

  • make_umi_qc_tables.sh

Interpretations

Duplex families are valuable for their low noise rate after collapsing, thus we'd like to see as high of a duplex "saturation" as possible. If this value is lower, we may not have captured enough of the original molecules to find both strands after PCR replication.

Theoretical Method

In this plot we investigate the number of families of each discrete size for simplex reads, which consist of 3 or more read pairs from one of the two strands.

Technical Methods

• Tools Used:

Theoretical Method

Coverage of each genomic region in the ACCESS panel is grouped on a per-sample basis, and a distribution of these values is plotted. Each sample is normalized by the median coverage value of that sample to align all peaks with one another and correct for sample-level differences.

Technical Methods

• Tool Used:

Theoretical Method

Count the number of positions in the bam that have an alt allele frequency of >0 and <2%

Theoretical Method

Similarly to the Pool A Targets, coverage is calculated for each UMI family type, over the Pool B genomic bait regions. These coverage values are lower for cfDNA samples (which use a 50:1 pool ratio) but should be comparable for buffy coat samples (which use a 1:1 pool ratio).

Technical Methods

• Tool Used:

Theoretical Method

Noise is calculated in the following manner:

Our current threshold for this calculation is set to 2%. Therefore it should be noted that there may be certain noisy positions which are wrongfully excluded, and other sites with low-level true mutations which are wrongfully included in the calculation.

In addition, inserted bases will be included in this calculation, but neither deletions, nor masked bases (N) are considered as alt alleles, nor are they counted towards the total depth.

Technical Methods

• Tool Used:

  • Marianas

  • Waltz PileupMetrics

Interpretations

Noise level can be influenced by a number of factors, including sequencing depth (and therefore coverage), duplex family saturation, and tumor content. We normally see the noise level for Duplex bams in the Pool A regions to be less than .001% (when using a 2% threshold for positions that should be included in the calculation). This threshold is indicated by the yellow dotted line in the graph. Noise higher than this value might be an indicator of a sample processing issue.

Theoretical Method

Major contamination plot is a bar plot of the fraction of heterozygous positions per sample and is done to see if a patient’s sample is contaminated with DNA from an unrelated individual. This analysis also done using the ‘fingerprint’ SNPs in the panel. A SNP is considered heterozygous if the minor allele fraction is > 0.1.

The fraction of heterozygous positions in the sample is found using the formula below:

Technical Methods

• Tool Used:

  • Waltz PileupMetrics

  • fingerprinting.py

Interpretations

The fraction of heterozygous positions should be around 0.5. If the fraction is greater than 0.6, it is is considered to have major contamination.

Theoretical Method

For each position that crosses the noise threshold (usually set at 2%), base changes are counted for each of the 6 possible substitution types.

Technical Methods

• Tool Used:

  • Marianas

  • Waltz PileupMetrics

Interpretations

ACCESS cfDNA samples usually exhibit larger noise values for C>T transitions, possibly due to cytosine deamination. However, differences between samples are not unexpected. Our threshold for ACCESS samples is 0.001 (past which we would fail a sample).

Theoretical Method

The sample mix-up heatmap is used to identify any potential mispaired samples within the run. The analysis makes use of the >300 ‘fingerprint’ single nucleotide polymorphisms (SNPs) that are distributed throughout the genome. These SNPs included the 31 SNPs that are in Target Pool A and >250 SNPs located in the tiling probes in Target Pool B. Pairwise comparisons of these SNP sites are done against all samples in the run. Sites, where both samples are homozygous are identified and percent discordance is calculated using the formula below:

where homozygous mismatches are sites that are homozygous in both Reference and Query but do not match each other.

If there are <10 common homozygous sites, the discordance rate can not be calculated since this is a strong indication that coverage is too low and the samples failed other QC.

Any samples with a discordance rate of 5% or higher are considered mismatches.

Technical Methods

• Tool Used:

  • Waltz PileupMetrics

  • fingerprinting.py

Interpretations

Dark blue indicates a match. Samples from the same patients are expected to match.

Theoretical Method

Similarly for the Simplex read pairs, we investigate the number of families of each discrete size for duplex reads, which consist of fragments with at least 1 read pair mapping on each of the top and bottom strands.

Technical Methods

• Tools Used:

  • Marianas

  • make_umi_qc_tables.sh

• Input

  • collapsed_R1_.fastq

  • collapsed_R2_.fastq

• Output

  • family-sizes.txt

Interpretations

We expect duplex family size peak between 5 and 15 read pairs, which gives us confidence that there are enough unique molecules for adequate error correction during the collapsing process.

Theoretical Method

Minor contamination check is done to see if a patient’s sample is contaminated with little DNA from another unrelated individual. This analysis is done using the ‘fingerprint’ SNPs identified in the .

Theoretical Method

Expected Matches are extracted from the title file provided for the pipeline run. Any samples with the same Patient ID in the title file are expected to match. The Expected matches that were extracted from the title file are printed in the ExpectedMatches.txt file in the QC_Results/FPResults folder from the pipeline.

The pairs of samples are assigned their “Status” based on the following conditions:

• Expected Match: Expected to match from Title file and discordance rate<5% .

Waltz Metrics Files

• .covered-regions

  • chr

  • start

  • end

  • length

  • average coverage in the contiguous region

  • total coverage in the contiguous region

• .read-counts

  • bam file name

  • total reads

• .fragment-sizes

  • fragment-size

  • total frequency

• -pileup-without-duplicates.txt

  • similar to above but only unique fragments are counted

• -intervals.txt Header

  • chr

  • start

• -intervals-without-duplicates.txt

  • similar to above but only unique fragments are considered

After aggregate_bam_metrics.sh (aggregate across samples):

• waltz-coverage.txt - per sample coverage calculated across chosen genomic intervals

• fragment-sizes.txt - fragment size distributions for all samples

Theoretical Method

Minor contamination from duplex bams for tumor sample (identified by the title file) is additionally checked. This analysis is done using the same fingerprint SNP in identified in the FP_configuration file although there is a 200x coverage threshold.

The Minor Contamination Rate is the average (mean) minor allele frequency from homozygous fingerprint SNPs, where homozygous sites as those harboring < 5% minor allele frequency in the sequence data.

Technical Methods

• Tool Used:

  • Waltz PileupMetrics

  • fingerprinting.py

Interpretations

Samples with Duplex Minor contamination rates of >0.002 are considered contamination.

Theoretical Method

Extract read counts for mutation hotspots from "normal" sample pileups. Then look into tumor samples to determine whether these mutations may have been due to contamination of tumor into the normal. Unfiltered bams are used for the normal samples to widen the search for hotspots, and duplex bams are then used for tumor samples.

Technical Methods

• Tool Used:

  • Waltz PileupMetrics

  • BioinfoUtils.jar

  • plots_module.r

• Input

  • sample_id-duplex-pileup.txt (for duplex noise calculation)

• Output

  • hotspots-in-normals.txt

Interpretations

In the provided example we can see that there was potential contamination of the tumor sample into the normal sample for C-P835W4, as indicated by the 7 unfiltered reads that matched a mutation from the tumor. This may be due to improper separation of tumor and normal sample during extraction, or clonal hematopoiesis.

Theoretical Method

Sex is inferred by looking at the average coverage for Tiling_SRY_Y:2655301 and Tiling_USP9Y_Y:14891501 probes in the All Unique bams (found from the intervals file in the Waltz output for Pool B). When the sum of the average coverage per interval (2 on Y) is greater that 50, the sample is classified as male. If the inferred sex does not match the reported sex, it is classified as a mismatch. Reported sex is from the title file.

Technical Methods

• Tool Used:

  • Waltz PileupMetrics

  • fingerprinting.py

Interpretations

Sex mismatches are an indication of a sample mixup. Low coverage, especially in the Y Chromosome may lead to a false positive.