Skip to content

QC a sequencing run

Goal: turn a raw sequencing sample and a reference genome into a single, quantitative quality-control (QC) report.

snipe qc is alignment-free: it scores a sample against a reference from k-mer sketches, with no aligner and no BAM. The metrics it produces — coverage, per-position depth, error and mutation rates, contamination, and sex-chromosome signals — cover the same questions an alignment-based QC tool such as qualimap answers from a BAM, computed here from a .snipesig sketch instead.

You will learn:

  • how to sketch a reference (--genome) and a sample (--sample) into .snipesig signatures;
  • how to run snipe qc to score the sample against the reference;
  • the three metric visibility tiers — default, --advanced, and --hidden;
  • how to read the headline metrics: coverage, sequencing-error rate, and mutation rate.
  1. Sketch the reference genome.

    A reference is sketched with --genome, which stores per-chromosome coordinates alongside the edgemers. It accepts exactly one FASTA:

    Terminal window
    snipe sketch --genome reference.fasta -o reference.snipesig

    This uses the defaults --k1-size 51, --k2-size 53, and --scale 10000. The sample you QC against it must be sketched with the same k1, k2, and scale.

  2. Sketch the sample.

    The reads (or an assembled sample) are sketched with --sample:

    Terminal window
    snipe sketch --sample sample.fastq.gz -o sample.snipesig

    --sample can be repeated to fold several files into one signature. Inputs may be .gz-compressed FASTA/FASTQ.

  3. Run the QC.

    Point snipe qc at the reference and the sample, and choose an output TSV. Both --reference and -o/--output are required:

    Terminal window
    snipe qc --reference reference.snipesig --sample sample.snipesig -o qc.tsv

    --sample may be repeated (or use --samples-list to read paths from a file) to score many samples in one run — each becomes a row in qc.tsv.

  4. Read qc.tsv.

    The output is tab-separated with dozens of columns — one metric per column, one sample per row. Here is a real run of the phiX174 genome scored against itself (its header plus first data row):

    qc.tsv (first row)
    experiment_id	k-size	edgemer_extension_length	sketching_scale	sketch_file_path	no_of_input_sequences	avg_sequence_length	GC%	no_of_input_bases	no_of_input_kmers	kmer_yield	fracminhash_precision	basepair_change_rate	sequencing_error_rate	mutation_rate	total_reference-mapped_bases	mean_depth_of_reference_coverage	reference_mapping_rate	fraction_of_reference_covered_by_sample_(%)	mean_depth_at_1x_covered_reference_bases	total_ROI-mapped_bases	mean_depth_of_ROI_coverage	ROI_mapping_rate	fraction_of_ROI_covered_by_sample_(%)	mean_depth_at_1x_covered_ROI_bases	ROI_enrichment_score	chromosomal_mean_abundance_CV	chrX_ploidy	chrY_coverage	predicted_contamination_index
    phix	51	2	1	phix.snipesig	1	5386.0	44.7642	5386	5336	0.9907	1	0.000000	0.000000	0.000000	5386	1.0094	1.000000	100	1.0094	0	0	0.000000	0	0	0	0	0	0	0
    

    Because this sample is the reference, the ideal-case numbers are what you’d expect: reference_mapping_rate = 1.000000, fraction_of_reference_covered_by_sample_(%) = 100, mean_depth_of_reference_coverage1.0094, and both sequencing_error_rate and mutation_rate = 0.000000. On real data these move away from the ideal, and that movement is the signal.

snipe qc groups every metric into a visibility tier. Add a flag to widen the report:

Tier Flag What it adds
Default (none) The everyday metrics — always written to the TSV.
Verbose --advanced Extra diagnostic columns for deeper inspection.
Hidden --hidden Developer-only metrics (intentionally kept out of qc --help).

The flags are additive — the columns above are the default tier, and each flag appends more columns to every row.

  • Coveragemean_depth_of_reference_coverage is the average number of times each covered reference position was seen, and fraction_of_reference_covered_by_sample_(%) (a fraction in [0, 1] via the API, a percentage in the TSV) tells you how much of the reference the sample touched at all. Low breadth with high depth means an uneven or targeted run.
  • sequencing_error_rate — estimated from edgemers whose K1 matched the reference but whose K2 extension mismatched and occurred only once (singletons). Single-occurrence changes are more likely miscalled bases than real variants, so a higher value indicates a noisier run.
  • mutation_rate — estimated from edgemers whose K2 extension mismatched and recurred (non-singletons). Multi-occurrence changes are more likely true biological variants than sequencing errors.
  • basepair_change_rate — the combined per-base change rate that the error and mutation rates decompose, corrected for the K1-to-K2 extension length.

For the full catalog — every column, its definition, and its tier — see the QC metrics reference, and read QC metrics explained for the intuition behind the error-vs-mutation split.

The same workflow through the Python API:

import snipe
reference = snipe.sketch_reference("reference.fasta")
sample = snipe.sketch(["sample.fastq.gz"])
metrics = sample.qc(reference=reference) # pandas Series, indexed by metric name
print(metrics["reference_mapping_rate"])
print(metrics["sequencing_error_rate"])

Signature.qc() returns a pandas Series indexed by metric name; pass advanced=True or hidden=True to widen it to match the CLI tiers. Metrics whose names end in _(%) come back as fractions in [0, 1]. You can also load a signature you already sketched with snipe.load("sample.snipesig").

You sketched a reference with --genome and a sample with --sample at matching k1/k2/scale, scored the sample with snipe qc --reference ... --sample ... -o qc.tsv, learned the default/--advanced/--hidden tiers, and read coverage, error rate, and mutation rate from the tab-separated result.