fastcat

fastcat summarises FASTQ and BAM data and can create offline quality-control reports. Use fastcat qc for new analyses; it supports read filtering, optional barcode demultiplexing, alignment metrics, and coverage summaries where applicable.

fastcat is also quite fast! See Performance.

Installation

Conda

A conda package is available through our conda channel:

mamba create -n fastcat -c conda-forge -c bioconda -c nanoporetech fastcat

Build from source

fastcat depends on bundled third-party sources in git submodules (third_party/htslib, third_party/zlib-ng, and third_party/argp-standalone). Purrflate is an optional acceleration submodule for multithreaded fastq.gz decompression.

Building requires several prerequisites:

• git (with submodule support)
• C toolchain: gcc (or clang), make, bash, patch
• Autotools: autoconf, automake (provides autoreconf/autoheader)
• Development libraries: zlib, bzip2, xz/lzma, libcurl, openssl
• Cmake is required for the optional depenency of libdeflate

For example on Debian or Ubuntu systems run:

sudo apt update
sudo apt install -y git gcc make bash patch autoconf automake \
  zlib1g-dev libbz2-dev liblzma-dev libcurl4-gnutls-dev libssl-dev

Then to clone and build fastcat run:

git clone --recursive https://github.com/epi2me-labs/fastcat.git
cd fastcat
make fastcat

--recursive initialises the public dependencies and deliberately skips the private optional Purrflate submodule. If you already cloned without submodules, run the following to initialise the public dependencies:

git submodule update --init third_party/htslib third_party/zlib-ng third_party/argp-standalone

Internal builds that have the Purrflate submodule use it automatically. Set NOPURRFLATE=1 to explicitly select the portable zlib FASTQ reader, while retaining Fastcat's normal threading support.

CLI overview

The command-line interface is split into multiple sub-commands with distinct tasks:

$ ./fastcat --help
Usage: ./fastcat <subcommand> [options]
Subcommands:
  qc           Unified FASTQ/BAM quality-control and reporting mode
  fastq        (Deprecated) FASTQ concatenation/statistics/filtering mode
  bam          (Deprecated) BAM alignment statistics mode
  lint         Standalone FASTQ low-complexity filter
  index        BAM chunk index tools (bamindex)
  json-schema  Write the summary.json JSON Schema to standard output
Try './fastcat <subcommand> --help' for subcommand options.

fastcat qc

Use fastcat qc to summarise FASTQ or BAM data and, optionally, create an offline HTML report. fastcat detects a common FASTQ or BAM format from named input files, directories, and file-of-filenames. Specify --input-fmt when reading input data from standard input, or to require a particular format.

# FASTQ quality-control report
fastcat qc reads.fastq.gz --html-report --output fastcat-qc

# BAM quality-control report
fastcat qc alignments.bam --html-report --output fastcat-qc

The primary QC output of fastcat qc is a summary.json file which includes all information garnered from the input files.

FASTQ and unaligned BAM inputs can be filtered with the read-filtering options and are written to standard output. Use --reads-output PREFIX to save them as PREFIX.fastq.gz, or PREFIX.bam with --bam-out, beneath the report directory. When read filtering is enabled, read-level QC metrics describe reads that pass the filters; the summary also reports counts of records rejected by each filter. For a BAM with reference sequences, fastcat collects alignment statistics only; it does not filter or write reads. The --region and --bed options can be used to select genomic regions from indexed BAMs. With one coordinate-sorted reference BAM, fastcat includes reference coverage summaries. Use --coverage to additionally write position-level coverage BED.gz files. For BAM input, --poly-a enables poly-A metrics for alignments with suitable data.

Useful options include --html-report, --output DIRECTORY, --threads THREADS, --sample NAME, --demultiplex, and --legacy-outputs. Run fastcat qc --help for the complete option reference. Use --jsonl-progress PATH when a wrapper needs machine-readable lifecycle and progress events; each line in that file is a complete JSON object. Progress events are emitted at most once per second, independently of the less-frequent human-readable terminal progress log. Every event has schema_version, a per-run monotonically increasing sequence, an ISO 8601 UTC timestamp, and an event type; remaining fields are specific to that event.

Legacy commands

fastcat fastq and fastcat bam are deprecated compatibility commands. New QC analyses should use fastcat qc; the older commands remain available for read concatenation, filtering, and existing automation. They retain their historical text-output default, which can be disabled with --no-legacy-outputs. Use fastcat fastq --help or fastcat bam --help for their option reference.

Outputs

fastcat qc writes its results under -o/--output DIRECTORY (default fastcat-qc; the directory must not already exist). Results are first written to DIRECTORY.tmp and published as DIRECTORY only after a successful run. Use --keep-failed-output to retain that temporary directory for debugging after a failure. Input and output failures use the conventional sysexits status values where applicable; an interrupted run exits with 128 + its signal number after removing its temporary output. With --demultiplex, it writes one output directory for each observed barcode under <output>/{unclassified,barcodeXXXX/}.

Each output directory contains a summary.json file, a machine-readable summary of the analysis. The --html-report option adds an offline interactive report. With --legacy-outputs enabled, fastcat also writes legacy aggregate text outputs: per_file.tsv, runids.tsv, basecallers.tsv, histogram files, BAM flagstats.tsv, and coverage summary and distribution tables where applicable. Their aggregate information is included in summary.json. per_read.tsv, requested separately with --per-read, is a detailed export and is not controlled by --legacy-outputs. Enabling the --coverage option results in detailed coverage BED.gz traces being written, in the case of a single coordinate-sorted reference BAM.

Read output

For FASTQ and unaligned BAM input, fastcat qc writes filtered reads to standard output in FASTQ format, or unaligned BAM with --bam-out. It does not write a read stream for reference-bearing BAM input. With --demultiplex, reads are written to barcode-specific output directories instead.

summary.json

The JSON summary contains statistics about the input data, including read and base yield, read length, mean quality, filtering outcomes, sequencing-over-time metrics, channel metrics, and quality diagnostics. Reference-bearing BAM input also adds alignment and coverage summaries. Most values are directly readable JSON; larger arrays are encoded compactly for efficient storage.

The following documentation is generated automatically from a JSON schema document describing the file contents. The JSON Schema document can be obtained by running fastcat json-schema from the CLI. Use this file if you need a machine readable specification of the summary.json files.

Dataframe encoding

For large array-like data the summary.json files use a self-describing base64-encoding. The descriptions contain the data type, the shape of the data and finally the encoded data.

The table below enumerates various definitions seen in the JSON schema document describing summary.json files. Examples of these structures are given folowing the table.

For example, 1D and 2D arrays are represented as array objects:

{
  "one_dimensional": {
    "name": "my_1d_array", (optional)
    "kind": "array",
    "dtype": "<u2",
    "shape": [3],
    "data": "AQACAAMA"
  }
  "two_dimensional": {
    "kind": "array",
    "dtype": "<f4",
    "shape": [2, 2],
    "data": "AACAPwAAAEAAAEBAAACAQA=="
  }
}

To optimize filesizes, the data type size, as represented by the dtype field in the above examples, may differ across files for the same attribute. Parsers should therefore not assume that a given array or data frame column is always encoded with the same data type, but rather inspect the dtype field.

Tabular, dataframe objects, data are simply lists of named array objects:

{
  table: {
    "kind": "dataframe",
    "shape": [2, 2],
    "columns": [
      {
        "name": "chrom",
        "kind": "array",
        "dtype": "string",
        "shape": [2]
        "data": ["chr1", "chr2"]
      },
      {
        "name": "mean_depth",
        "kind": "array",
        "dtype": "<u4",
        "shape": [2]
        "data": "ZAAAAMgAAAA="
      }
    ]
  }
}

run_summary — Run summary

The principal human-readable overview for this output sample or barcode: accepted-read yield and count, filtering outcomes, sequencing-run metadata, and N/L read-length statistics.

sample — Sample metadata

Sample and barcode identifiers for this output.

files — Input file summaries

Per-input-file accepted-read summaries. In a demultiplexed run, each row is scoped to the barcode represented by this summary.json file.

counts — Aggregate count summaries

Machine-readable aggregate accepted-read and base totals, with frequency counters that preserve the number of reads associated with each observed run ID and basecaller model.

channels — Per-channel metrics

Sparse channel-and-mux aggregates used for the channel-occupancy maps: accepted read and base totals, plus mean quality and, for aligned input, mean alignment concordance. This value is null when usable channel metadata is unavailable or when an output contains multiple flowcells.

The data dataframe columns are:

distributions — Read statistic histograms

Histogram data for canonical reads: primary mapped and primary unmapped BAM records share the same read-length and mean-quality distributions.

The read_length.data dataframe columns are:

The quality.data dataframe columns are:

time — Sequencing metrics over time

Chronological per-bin read and base yield, cumulative totals, and read-statistic quantiles. These values underpin the sequencing-over-time plots. This value is null when usable timing metadata is unavailable.

The data dataframe columns are:

heatmaps — Quality diagnostic heatmaps

Two-dimensional count matrices relating mean read quality to read length, translocation speed, or alignment concordance. Axis metadata identifies the metric and its numeric range.

alignment — Alignment summary metrics

BAM-only record-level mapping counts, canonical-read mapping totals, per-reference flag counts, and distributions for alignment concordance, aligned query coverage, and poly-A length. This section is omitted for FASTQ summaries.

The distributions.accuracy.data and distributions.query_coverage.data dataframe columns are:

The distributions.poly_a_length.data dataframe columns are:

coverage — Reference coverage

Reference-depth tables for a single coordinate-sorted BAM with reference sequences. Target sets may cover the whole reference, fixed-length segments, or named BED intervals; the summary and cumulative-distribution dataframes directly mirror fastcat's legacy coverage tables. This section is omitted for FASTQ, unaligned BAM, multiple BAM inputs, and unsorted aligned BAM.

The targets[].summary.data dataframe columns are:

The targets[].distribution.data dataframe columns are:

Histogram outputs

These outputs are considered legacy outputs and may be removed in future versions. It is recommended to use the summary.json output instead.

The program creates a histograms/ directory (per barcode, when demultiplexed) containing:

• length.hist - read length histogram.
• quality.hist - read mean base-quality score histogram.
• accuracy.hist (bam only) - read alignment accuracy histogram.
• coverage.hist (bam only) - read alignment coverage histogram (the proportion of each read spanned by the alignment and not clipped).
• polya.hist (bam only, with --poly-a) - poly-A tail length histogram.

The format of the histogram files is a tab-separated file of sparse, ordered intervals [lower, upper):

lower    upper    count

The final bin may be unbounded, which is signified by a 0 entry for the upper bin edge.

Legacy per-file outputs

per_file.tsv is a legacy per-input processing diagnostic, rather than part of the supported aggregate JSON interface. It is a tab-separated file with columns:

filename        n_seqs  n_bases min_length      max_length      mean_quality    f_file_ok       f_stream_error  f_qual_missing  f_qual_truncated        f_unknown_error r_record_seen   r_record_ok     r_too_long r_too_short     r_low_quality   r_dust_masked
test/data/bc0.fastq.gz  10      4599    191     965     9.86    1       0       0       0       0       10      10      0       0       0       0
...

where the mean_quality column is the mean of the per-read mean_quality values. The f_* fields detail file-level parse errors. r_record_seen is records parsed/considered, while r_record_ok is records that passed filters and were used downstream.

The basecallers.tsv and runids.tsv files contain per-file counts of unique basecallers and MinKNOW run identifiers observed in the input files. Their structures are, respectively:

filename        basecaller      count
test/data/bc0.fastq.gz          10
...

and

filename        run_id  count
test/data/bc0.fastq.gz  5a21d8a6996146deceeaea3784244c52741cae93        10

Per-read output

The (optional, enabled with -p/--per-read) per_read.tsv is a tab-separated file with columns:

read_id filename        runid   read_length     mean_quality    channel read_number     start_time
32e13a1c-4171-4706-b6ce-a32c0f65fa16    test/data/bc0.fastq.gz  5a21d8a6996146deceeaea3784244c52741cae93        326     12.31   282     9       2021-04-20T17:00:40Z
b87f011e-b802-4993-8f56-fd240b2e784f    test/data/bc0.fastq.gz  5a21d8a6996146deceeaea3784244c52741cae93        407     9.13    213     19      2021-04-20T17:00:41Z
6f64aedb-bb8e-4777-b494-43e661841e06    test/data/bc0.fastq.gz  5a21d8a6996146deceeaea3784244c52741cae93        355     9.98    67      13      2021-04-20T17:00:41Z
c372fb2c-dd45-4feb-81b2-c167c3d1ce93    test/data/bc0.fastq.gz  5a21d8a6996146deceeaea3784244c52741cae93        317     8.14    337     18      2021-04-20T17:00:41Z
18d04e8d-2816-4986-8e1b-e5be676837fc    test/data/bc0.fastq.gz  5a21d8a6996146deceeaea3784244c52741cae93        965     7.57    507     18      2021-04-20T17:00:41Z
...

(bam only) additional alignment columns are appended: ref, coverage, ref_coverage, qstart, qend, rstart, rend, aligned_ref_len, direction, length, match, ins, del, sub, iden, acc, duplex.

The mean quality is defined as:

-10 * log10(mean(10^(Q/-10)))

where Q are the set of all per-base quality scores for the read. In many cases (where relevant metadata is found in the read header) the mean quality will not be recomputed but taken from the basecaller's value, which may be subtly different. Recalculation of the mean quality can be forced with --recalc-qual.

Alignment flag outputs — bam only

This output is considered a legacy output and may be removed in future versions. It is recommended to use the summary.json output instead.

The file flagstats.tsv provides similar output to samtools flagstats, enumerating counts for different SAM flags for each reference sequence:

ref     sample_name     total   primary secondary       supplementary   unmapped        qcfail  duplicate       duplex  duplex_forming
cneofor_1       0       0       0       0       0       0       0       0       0
cneofor_2       0       0       0       0       0       0       0       0       0
cneofor_3       1       0       0       1       0       0       0       0       0
ecoli1          373     347     5       21      0       0       0       0       0
ecoli2          0       0       0       0       0       0       0       0       0

Alignment coverage outputs — bam only

The JSON summary contains the contents of the {name}.summary.txt and {name}.dist.txt files referenced below. It does not contain the complete contents of the BED formatted files.

The program can optionally output coverage information (--coverage, provided the input BAM is coordinate sorted). The outputs are inspired by those from mosdepth but have been adapted to be more readily interpretable in line with user expectations. A set of files is produced for each user-provided BED file (--coverage-beds/--coverage-names), as well as for the whole genome. The files are:

• {name}.bed.gz - a bgzip-compressed BED file with position-level accuracy.
• {name}.fwd.bed.gz - as above but only for forward strand alignments.
• {name}.rev.bed.gz - as above but only for reverse strand alignments.
• {name}.summary.txt - a tab-separated summary of coverage for each region in the BED file. And a final total line for all regions in the BED.
• {name}.dist.txt - a tab-separated file listing the cumulative distribution function of coverage across all regions in the BED file (i.e. "what proportion of positions are covered at at least X depth?").

The summary file contains user-configurable entries (--thresholds) enumerating a sparse coverage CDF: the default is to output the proportion of positions covered at 1x, 5x, 10x, 20x, 30x, and 40x. This is done independently for all regions in the BED, as well as the final total line.

Performance

fastcat qc is somewhat highly optimised. It is designed as a single-pass streaming analysis: records are parsed, filtered, and passed to the requested metric sinks without retaining the read set in memory. The FASTQ reader uses our high-performance, zero-copy, parallel gzip decoder purrflate. BAM input and BAM output use HTSlib's parallel functionality.

When unaligned BAM output is requested from FASTQ, Fastcat uses a specialised converter rather than the general-purpose HTSlib record constructor. Fastcat also trusts a valid qs/QS mean-quality tag by default, avoiding a per-base quality scan.

For coordinate-sorted BAM coverage, Fastcat records interval-boundary deltas in separate forward- and reverse-strand arrays for one contig at a time. A cumulative sweep then produces coverage summaries and, when requested, coalesced bedGraph intervals.

Several distribution and time-series summaries use KLL quantile sketches. Their compactors retain a bounded O(k log n) sample of n observations instead of every observed value, allowing Fastcat to estimate quantiles and smooth report densities with predictable memory use even for very large runs.

Runtime numbers are shown below for a dataset containing 500,000 reads with an average length of 20 kbase. Results are presented for four input scenarios:

• BAM file (including methylation data), 12Gb
• FASTQ.gz file, 8.2Gb
• FASTQ.gz data piped from standard input stream
• FASTQ plain-text file, 18Gb

The program was run as:

  fastcat qc --input-fmt fastq --metrics-only --threads N --html-report --output REPORT

Production distributions of fastcat benefit from the use of our purrflate parallel gzip decoder for FASTQ input; building from source will not show scaling with thread count for FASTQ.

The plain-text FASTQ control shows baseline performance of the fastcat code. The results show that for compressed inputs (including BAM) the runtime is overwhelmingly dominated by decompression performance. The BAM results are slower than FASTQ as the files are larger, requiring decompression of the additional data.

fastcat lint

fastcat lint is a lightweight FASTQ filter for removing reads dominated by low-complexity sequence. It applies the SDUST algorithm to each read, writes reads at or below the configured masked-base proportion to standard output, and reports discarded reads on standard error. It is useful as a small pipeline stage when low-complexity reads need to be excluded before downstream processing.

$ ./fastcat lint --help
Usage: lint [OPTION...] <reads.fastq>
fastlint -- apply sdust algorithm to input files.

 General options:
  -p, --max-proportion=PROPORTION
                             Maximum allowable proportion of masked bases in a
                             read to keep the read (default: 0.95).
  -t, --threshold=THRESHOLD  Threshold for repetition (default: 20).
  -w, --window=WINDOW        Window size (default: 64).

  -?, --help                 Give this help list
      --usage                Give a short usage message
  -V, --version              Print program version

fastcat index

fastcat index splits an unaligned, unsorted BAM logically into record-count chunks without rewriting the source file. build creates a .bci index of BGZF offsets, fetch streams one indexed chunk as BAM to standard output, and dump makes the index contents inspectable as text. This is intended for scattering record-based BAM work across parallel workflow tasks; it is not a coordinate index and does not keep records for a query together.

$ ./fastcat index --help
* bamindex help          Print general help or help about a subcommand.
* bamindex build         Build a BAM index.
* bamindex fetch         Fetch records from a BAM using an index.
* bamindex dump          Dump an index fetch to text.

fastcat index build

$ ./fastcat index build --help
Usage: build [OPTION...] <reads.bam>
bamindex build -- create a BAM index corresponding to batches of records.

 General options:
  -c, --chunk-size=SIZE      Number of records in a chunk.
  -t, --threads=THREADS      Number of threads for BAM processing.

  -?, --help                 Give this help list
      --usage                Give a short usage message
  -V, --version              Print program version

fastcat index fetch

$ ./fastcat index fetch --help
Usage: fetch [OPTION...] <reads.bam.bci>
bamindex fetch -- fetch records from a BAM according to an index.

 General options:
  -c, --chunk=SIZE           Chunk index to retrieve.
  -t, --threads=THREADS      Number of threads for BAM processing.

  -?, --help                 Give this help list
      --usage                Give a short usage message
  -V, --version              Print program version

fastcat index dump

$ ./fastcat index dump --help
Usage: dump [OPTION...] <reads.bam.bci>
bamindex dump -- dump a BAM chunk index to stdout as text.

  -?, --help                 Give this help list
      --usage                Give a short usage message
  -V, --version              Print program version