[](https://github.com/mbhall88/rasusa/actions/workflows/rust-ci.yaml)
[](https://opensource.org/licenses/MIT)
[](https://github.com/mbhall88/rasusa/releases)
[](https://doi.org/10.21105/joss.03941)
[](https://doi.org/10.46471/gigabyte.180)
**Ra**ndomly **su**b**sa**mple sequencing reads or alignments.
> Furqon, A. D. C., Roberts, L. W., & Hall, M. B. (2026). Efficient downsampling of genome alignments with Rasusa. Gigabyte, 2026, gigabyte180. https://doi.org/10.46471/gigabyte.180
> Hall, M. B., (2022). Rasusa: Randomly subsample sequencing reads to a specified coverage. Journal of Open Source
> Software, 7(69), 3941, https://doi.org/10.21105/joss.03941
[TOC]: #
## Table of Contents
- [Table of Contents](#table-of-contents)
- [Install](#install)
- [Usage](#usage)
- [Benchmark](#benchmark)
- [Contributing](#contributing)
- [Citing](#citing)
## Install
### Precompiled binary
The quickest way to get `rasusa` is with the install script:
```shell
curl -sSL https://github.com/mbhall88/rasusa/releases/latest/download/install.sh | sh
# or with wget
wget -nv -O - https://github.com/mbhall88/rasusa/releases/latest/download/install.sh | sh
```
The installer will identify the architecture and OS, download the latest binary for your system, and move it to `/usr/local/bin` (on Linux/macOS).
If you would like to install it to a different location, you can pass the `-b` or `--bin-dir` option to the script.
```shell
curl -sSL https://github.com/mbhall88/rasusa/releases/latest/download/install.sh | sh -s -- -b /my/custom/bin
```
You can also pass options to the script like so
```
$ curl -sSL https://github.com/mbhall88/rasusa/releases/latest/download/install.sh | sh -s -- --help
install.sh [option]
Fetch and install the latest version of rasusa, if rasusa is already
installed it will be updated to the latest version.
Options
-V, --verbose
Enable verbose output for the installer
-f, -y, --force, --yes
Skip the confirmation prompt during installation
-p, --platform
Override the platform identified by the installer [default: apple-darwin]
-b, --bin-dir
Override the bin installation directory [default: /usr/local/bin]
-a, --arch
Override the architecture identified by the installer [default: x86_64]
-B, --base-url
Override the base URL used for downloading releases [default: https://github.com/mbhall88/rasusa/releases]
-h, --help
Display this help message
```
### `cargo`
[](https://crates.io/crates/rasusa)
```sh
cargo install rasusa
```
Alternatively, for more modern Rust-based binary management, we recommend using [`cargo-binstall`][binstall]:
```shell
cargo binstall rasusa
```
[binstall]: https://github.com/cargo-bins/cargo-binstall
### `conda`
[](https://anaconda.org/bioconda/rasusa)
[](https://anaconda.org/bioconda/rasusa)
Prerequisite: [`conda`][conda] (and bioconda channel [correctly set up][channels])
```sh
conda install rasusa
```
### Container
Docker images are hosted on [GitHub Container Registry][ghcr].
#### `apptainer`
Prerequisite: [`apptainer`][apptainer]
```sh
URI="docker://ghcr.io/mbhall88/rasusa"
apptainer exec "$URI" rasusa --help
```
The above will use the latest version. If you want to specify a version then use a
[tag][ghcr] (or commit) like so.
```sh
VERSION="2.2.2"
URI="docker://ghcr.io/mbhall88/rasusa:${VERSION}"
```
#### `docker`
Prerequisite: [`docker`][docker]
```sh
docker pull ghcr.io/mbhall88/rasusa
docker run ghcr.io/mbhall88/rasusa --help
```
You can find all the available tags on the [container registry][ghcr]. Note: versions
prior to 0.4.0 were housed on [Docker Hub](https://hub.docker.com/r/mbhall88/rasusa).
And versions from 0.4.0 to 2.2.0 were on [quay.io][quay.io].
### `homebrew`
Prerequisite: [`homebrew`][homebrew]
```sh
brew install rasusa
```
### Build locally
Prerequisite: [`rust` toolchain][rust]
```sh
git clone https://github.com/mbhall88/rasusa.git
cd rasusa
cargo build --release
target/release/rasusa --help
# if you want to check everything is working ok
cargo test --all
```
## Usage
### Basic usage - reads
Subsample FASTQ reads or unaligned SAM/BAM/CRAM.
```
rasusa reads --coverage 30 --genome-size 4.6mb in.fq
```
`rasusa reads` also supports unaligned SAM/BAM/CRAM input:
```
rasusa reads --coverage 30 --genome-size 4.6mb in.bam
```
The above commands will output the subsampled file to `stdout`.
Or, if you have paired Illumina
```
rasusa reads --coverage 30 --genome-size 4g -o out.r1.fq -o out.r2.fq r1.fq r2.fq
```
For more details on the above options, and additional options, see below.
### Basic usage - alignments
Subsample alignments
```
rasusa aln --coverage 30 in.bam | samtools sort -o out.bam
```
this will subsample each position in the alignment to 30x coverage.
### Required parameters
`rasusa` has three required options for the `reads` command, and two required options for the `aln` command.
#### Input
This positional argument specifies the file(s) containing the reads or unaligned alignments you would like to subsample.
For the `reads` command, the file(s) can be in FASTA, FASTQ, or unaligned SAM, BAM, or CRAM format. FASTA and FASTQ files can be compressed (with a tool such as `gzip`).
If two files are passed to `reads`, `rasusa` will assume they are paired-end reads.
For the `aln` command, the file must be a valid coordinate-sorted SAM/BAM/CRAM file.
> Bash wizard tip 🧙: Let globs do the work for you `r*.fq`
#### Coverage
##### `-c`, `--coverage`
> Not required if [`--bases`](#target-number-of-bases) is present for `reads`
This option is used to determine the minimum coverage to subsample the reads to. For the `reads` command, it can
be specified as an integer (100), a decimal/float (100.7), or either of the previous
suffixed with an 'x' (100x). For the `aln` command, it is an integer only.
Due to the method for determining how many bases are required to achieve the
desired coverage in the `reads` command, the actual coverage, in the end, could be slightly higher than
requested. For example, if the last included read is very long. The log messages should
inform you of the actual coverage in the end.
For the `aln` command, the coverage is the maximum number of reads that should be present at each position in the
alignment. If a position has fewer than the requested number of reads, all reads at that position will be included.
> [!NOTE]
> **For paired-end data:**
> To ensure 100% pair retention (i.e., no orphan reads) during subsampling, the `aln` command uses a two-pass strategy.
> It first subsamples the first segment (Read 1) at half the target coverage (`coverage / 2`) and then recovers
> the corresponding mates (Read 2).
>
> Because this approach is constrained by pairing rather than per-position depth alone, the resulting coverage
> is not guaranteed to be perfectly uniform across all positions. Instead, some fluctuations around the requested
> depth should be expected.
>
> If the requested coverage is a strict *minimum* requirement, we recommend setting `--coverage` slightly higher to account
> for these fluctuations.
> See the discussion in [this PR](https://github.com/mbhall88/rasusa/pull/118) for additional details.
#### Genome size
##### `-g`, `--genome-size`
> Not valid for `aln`
> Not required if [`--bases`](#target-number-of-bases) is present for `reads`
The genome size of the input is also required. It is used to determine how many bases
are necessary to achieve the desired coverage. This can, of course, be as precise or
rough as you like.
Genome size can be passed in many ways. As a plain old integer (1600), or with a metric
suffix (1.6kb). All metric suffixes can have an optional 'b' suffix and be lower, upper,
or mixed case. So 'Kb', 'kb' and 'k' would all be inferred as 'kilo'. Valid metric
suffixes include:
- Base (b) - multiplies by 1
- Kilo (k) - multiplies by 1,000
- Mega (m) - multiplies by 1,000,000
- Giga (g) - multiplies by 1,000,000,000
- Tera (t) - multiplies by 1,000,000,000,000
Alternatively, a [FASTA/Q index file][faidx] can be given and the genome size will be
set to the sum of all reference sequences in it.
> [!TIP]
> If you want to use `rasusa` in a scenario where you don't know what the genome size is,
> such as in an automated pipeline that can take in any kind of organism, you could estimate
> the genome size with something like [`lrge`](https://github.com/mbhall88/lrge) (#shamelessplug).
>
> ```
> $ gsize=$(lrge reads.fq)
> $ rasusa reads -g $gsize -c 10 reads.fq
> ```
> `lrge` is designed for long reads. If you want to estimate the genome size from short
> reads, you could use something like [Mash](https://github.com/marbl/Mash) or
> [GenomeScope2](https://github.com/tbenavi1/genomescope2.0). See [the `lrge` docs](https://github.com/mbhall88/lrge?tab=readme-ov-file#alternatives)
> for examples of how Mash/GenomeScope2 can be used for this task.
[faidx]: https://www.htslib.org/doc/faidx.html
### Optional parameters
#### Output
##### `-o`, `--output`
**`reads`**
> [!IMPORTANT]
> This parameter is required if passing paired Illumina data to `reads`.
By default, `rasusa` will output the subsampled file to `stdout` (if one file is given).
If you would prefer to specify an output file path, then use this option.
Output for Illumina paired files must be specified using `--output` twice - `-o out.r1.fq -o out.r2.fq`
The ordering of the output files is assumed to be the same as the input.
`rasusa reads` will attempt to automatically infer the output format and whether compression of the output
file(s) is required. It does this by detecting any of the supported extensions:
- `.fa`/`.fasta`: FASTA format
- `.fq`/`.fastq`: FASTQ format
- `.bam`: BAM format
- `.cram`: CRAM format
- `.sam`: SAM format
- `.gz`: will compress the output with [`gzip`][gzip]
- `.bz` or `.bz2`: will compress the output with [`bzip2`][bzip]
- `.lzma`: will compress the output with the [`xz`][xz] LZMA algorithm
- `.zst`: will compress the output with [`zstd`][zstd]
> [!NOTE]
> If no extension is matched, the output will be in the same format as the input. You can convert between formats
> by specifying the desired extension (e.g., input BAM and output FASTQ). However, you cannot convert from FASTA/FASTQ
> to unaligned SAM/BAM/CRAM--though you can convert from unaligned SAM/BAM/CRAM to FASTA/FASTQ.
**`aln`**
For the `aln` command, the output file format will be the same as the input if writing to stdout, otherwise it will be
inferred from the file extension.
> [!NOTE]
> The output alignment will most likely **not be sorted**. You can use `samtools sort` to sort the output. e.g.,
>
> ```
> rasusa aln -c 5 in.bam | samtools sort -o out.bam
> ```
[gzip]: http://www.gzip.org/
[bzip]: https://sourceware.org/bzip2/
[xz]: https://tukaani.org/xz/
#### Output compression/format
**`reads`**
`rasusa reads` provides two options for explicitly setting the output format and compression.
##### `-Z`, `--compress-type`
Use this option to manually set the compression algorithm to use for the output file(s).
It will override any format automatically detected from the output path. Note: this is only used for FASTA/FASTQ output.
Valid options are:
- `g`: [`gzip`][gzip]
- `b`: [`bzip2`][bzip]
- `l` or `x`: [`xz`][xz] LZMA algorithm
- `z`: [`zstd`][zstd]
- `u`: no compression
##### `-O`, `--output-format`
Use this option to manually set the output format to use for the output file(s).
It will override any format automatically detected from the output path.
Valid options are:
- `f` or `fasta`: FASTA
- `q` or `fastq`: FASTQ
- `s` or `sam`: SAM
- `b` or `bam`: BAM
- `c` or `cram`: CRAM
**`aln`**
##### `-O`, `--output-format`
Use this option to manually set the output file format for `rasusa aln`. By default, the same format as the input will be used, or the
format will be guessed from the `--output` path extension if given.
Valid options are:
- `s` or `sam`: SAM
- `b` or `bam`: BAM
- `c` or `cram`: CRAM
All values to this option are case insensitive.
#### Compresion level
##### `-l`, `--compress-level`
> `reads` only
Compression level to use if compressing the output. By default this is set to the default for the compression type being
output.
#### Target number of bases
##### `-b`, `--bases`
> `reads` only
Explicitly set the number of bases required in the subsample. This option takes the
number in the same format as [genome size](#genome-size).
> [!NOTE]
> If this option is given, genome size and coverage are not required.
#### Number of reads
##### `-n`, `--num`
> `reads` only
Explicitly set the number of reads in the subsample. This option takes the number in
the same format as [genome size](#genome-size).
When providing paired reads as input, this option will sample this many total read
pairs. For example, when passing `-n 20 r1.fq r2.fq`, the two output files will have
20 reads each, and the read ids will be the same in both.
*Note: if this option is given, genome size and coverage are not required.*
#### Fraction of reads
##### `-f`, `--frac`
> `reads` only
Explicitly set the fraction of total reads in the subsample. The value given to this
option can be a float or a percentage - i.e., `-f 0.5` and `-f 50` will both take half
of the reads.
> [!NOTE]
> If this option is given, genome size and coverage are not required.
#### One-pass (streaming) fraction subsampling
##### `-1`, `--one-pass`
> `reads` only
By default, `--frac` subsampling reads the input twice: once to measure it, then once
to subsample to an exact fraction. `--one-pass` instead reads the input exactly once,
deciding whether to keep each read (or, for paired/segmented input, each read
pair/template) independently with probability `--frac`, as it streams through. This is
faster and skips the first pass entirely, but the number of reads kept is only
approximately `--frac * total_reads`, not exact - see [`benches/README.md`][one-pass-accuracy]
for how close it lands in practice (close enough to ignore for a real sequencing run;
worth avoiding if you're subsampling a handful of reads).
> [!NOTE]
> `--one-pass` requires `--frac` - it isn't available for `--num`/`--bases`/`--coverage`,
> which need an exact read count or the input's total base count up front. It also
> cannot be combined with `--strict`, since the result is approximate by design.
For paired FASTA/Q, mates are always kept or dropped together. For a single unaligned
SAM/BAM/CRAM file, all of a segmented read's records are kept or dropped together as
one template - this requires the input to be grouped by read name (`samtools collate`,
not `samtools sort`). Name-sorted/grouped input is detected from the header
(`SO:queryname` or `GO:query`) where present; otherwise the first 50 records are
scanned for evidence of grouping, and the input is rejected up front if that scan finds
evidence it isn't.
#### Probability shorthand
##### `-p`, `--probability`
> `reads` only
Shorthand for `--frac <FLOAT> --one-pass` - `-p 0.1` is identical to `-f 0.1 -1`,
including producing byte-identical output for the same `--seed`. Values `>1` and `<=100`
are converted the same way as `--frac` - e.g., `-p 25` means `-p 0.25`.
#### Random seed
##### `-s`, `--seed`
This option allows you to specify the [random seed][seed] used by the random subsampler. By explicitly setting this
parameter, you make the subsample for the input reproducible. You only need to pass this parameter if you are likely
to want to subsample the same input file again in the future and want the same subset of reads. However, if you forget
to use this option, the seed generated by the system will be printed to the log output, allowing you to use it in the
future.
> [!WARNING]
> **v5.0.0 breaking change**: for `--num`/`--frac` subsampling, `rasusa` now selects the `k` reads to
> keep directly instead of shuffling every read index, which cuts selection time and memory from
> `O(n)` to `O(k)` (n = total reads, k = reads kept). This changes which reads a given `--seed`
> selects for `--num`/`--frac`. `--bases`/`--coverage` subsampling is unaffected - the same seed
> still produces identical output as in v4.x.
#### Strictly adhere to requested coverage
##### `-e`, `--strict`
> `reads` only
If the requested coverage, total bases, number of reads, or fraction of reads cannot be met, an error will be thrown.
By default, a warning is displayed, and the maximum possible coverage, total bases, number of reads, or fraction of reads is used.
#### Subsampling strategy
##### `--strategy`
> `aln` only
By default, `rasusa aln` uses the `stream` strategy, which implements a fast sweep-line algorithm with random priority.
It processes a coordinate-sorted alignment file in a single pass (two passes if using paired-end data) while maintaining an active set of reads in a heap,
ensuring that no position exceeds the target depth **N**.
This strategy provides the option `--swap-distance` (default: 5 bp), which limits the allowed distance when swapping
between reads encountered in the current scan and reads already in the heap.
Alternatively, users can select the `fetch` strategy. This approach repeatedly fetches overlapping reads,
shuffles them, and samples to the target depth **N**.
The fetch strategy provides additional controls:
- `--batch-size` (default: 10 kb): size of genomic window cached in memory
- `--step-size` (default: 100 bp): step size used when scanning along the chromosome to find overlapping reads
In most cases, the default `stream` strategy is recommended due to its speed and low memory usage.
See [this PR](https://github.com/mbhall88/rasusa/pull/118) for a discussion of performance and behavior differences between these two strategies.
#### Verbosity
##### `-v`
> `reads` only
Adding this optional flag will make the logging more verbose. By default, logging will
produce messages considered "info" or above (see [here][log-lvl] for more details). If
verbosity is switched on, you will additionally get "debug" level logging messages.
### Full usage
```text
$ rasusa --help
Randomly subsample reads or alignments
Usage: rasusa [OPTIONS] <COMMAND>
Commands:
reads Randomly subsample reads (FASTA/Q, unaligned SAM/BAM/CRAM)
aln Randomly subsample alignments to a specified depth of coverage
cite Get a bibtex formatted citation for this package
help Print this message or the help of the given subcommand(s)
Options:
-v Switch on verbosity
-h, --help Print help
-V, --version Print version
```
#### `reads` command
```text
$ rasusa reads --help
Randomly subsample reads (FASTA/Q, unaligned SAM/BAM/CRAM)
Usage: rasusa reads [OPTIONS] <FILE(S)>...
Arguments:
<FILE(S)>...
The FASTA/FASTQ or unaligned SAM/BAM/CRAM file(s) to subsample.
For paired Illumina, the order matters. i.e., R1 then R2. Single-file paired-end is also supported for unaligned SAM/BAM/CRAM.
Options:
-o, --output <OUTPUT>
Output filepath(s); stdout if not present.
For paired Illumina pass this flag twice `-o o1.fq -o o2.fq`
NOTE: The order of the pairs is assumed to be the same as the input - e.g., R1 then R2.
This option is required for paired input.
-g, --genome-size <size|faidx>
Genome size to calculate coverage with respect to. e.g., 4.3kb, 7Tb, 9000, 4.1MB
Alternatively, a FASTA/Q index file can be provided and the genome size will be set to the sum of all reference sequences.
If --bases is not provided, this option and --coverage are required
-c, --coverage <FLOAT>
The desired depth of coverage to subsample the reads to
If --bases is not provided, this option and --genome-size are required
-b, --bases <bases>
Explicitly set the number of bases required e.g., 4.3kb, 7Tb, 9000, 4.1MB
If this option is given, --coverage and --genome-size are ignored
-n, --num <INT>
Subsample to a specific number of reads
If paired-end reads are passed, this is the number of (matched) reads from EACH file. This option accepts the same format as genome size - e.g., 1k will take 1000 reads
-f, --frac <FLOAT>
Subsample to a fraction of the reads - e.g., 0.5 samples half the reads
Values >1 and <=100 will be automatically converted - e.g., 25 => 0.25
-1, --one-pass
Read the input exactly once, keeping each read independently with probability --frac, instead of measuring the input first for an exact result
Only supported for a --frac target: the number of reads kept is approximate (binomially distributed around --frac), not exact, so it cannot be combined with --strict, and it is not available for --num/--bases/--coverage targets, which need an exact count or the input's total base count up front.
Supports FASTA/Q (single-end or paired - paired mates are always kept or dropped together) and a single unaligned SAM/BAM/CRAM file (two separate SAM/BAM/CRAM files are not supported). For SAM/BAM/CRAM, a segmented read's records are grouped into one template by comparing each record only to the one immediately before it, which is only correct if the input is grouped by read name - name-sorted/grouped input (SO:queryname or GO:query in the header) is accepted outright, and anything else is checked by scanning the first 50 records for evidence of grouping (rejected if that scan finds a name that reappears after a different one; not a guarantee beyond those 50 records - collate the input first, e.g. with `samtools collate`, if in doubt).
-p, --probability <FLOAT>
Keep each read independently with this probability, streaming the input once
Shorthand for `--frac <FLOAT> --one-pass`: both spellings behave identically, including producing byte-identical output for the same seed. Values >1 and <=100 will be automatically converted, the same as --frac - e.g., 25 => 0.25.
The result is approximate (binomially distributed around the requested probability), not an exact fraction of the input.
-e, --strict
Exit with an error if the requested coverage/bases/reads is not possible
-s, --seed <INT>
Random seed to use
-v
Switch on verbosity
-Z, --compress-type <u|b|g|l|x|z>
u: uncompressed; b: Bzip2; g: Gzip; l: Lzma; x: Xz (Lzma); z: Zstd
Rasusa will attempt to infer the output compression format automatically from the filename extension. This option is used to override that. If writing to stdout, the default is uncompressed. Note: this is only used for FASTA/FASTQ output.
-O, --output-format <OUTPUT_FORMAT>
Explicitly set the output format.
If not provided, Rasusa will attempt to infer the format from the filename extension.
[possible values: fasta, fastq, sam, bam, cram]
-l, --compress-level <1-21>
Compression level to use if compressing output. Uses the default level for the format if not specified
-@, --threads <INT>
Number of threads to use for BAM (de)compression
Only relevant for BAM input/output - SAM is uncompressed and CRAM has no multithreaded codec, so this is a no-op for those formats and for FASTA/FASTQ. `--threads 1` is identical to not passing this option. Only BAM *reading* benefits here; writing BAM output remains single-threaded.
[default: 1]
-h, --help
Print help (see a summary with '-h')
-V, --version
Print version
```
#### `aln` command
```text
$ rasusa aln --help
Randomly subsample alignments to a specified depth of coverage
Usage: rasusa aln [OPTIONS] --coverage <INT> <FILE>
Arguments:
<FILE>
Path to the input alignment file (SAM/BAM/CRAM) to subsample
Note: An index (.bai) is required when using '--strategy fetch'.
Options:
-o, --output <FILE>
Path to the output subsampled alignment file. Defaults to stdout (same format as input)
The output is not guaranteed to be sorted. We recommend piping the output to `samtools sort`
-O, --output-format <FMT>
Output format. Rasusa will attempt to infer the format from the output file extension if not provided
-c, --coverage <INT>
The desired depth of coverage to subsample the alignment to
-s, --seed <INT>
Random seed to use
--strategy <STRATEGY>
Subsampling strategy
Possible values:
- stream: A linear scan approach using sweep line algorithm with random priority. Requires sorted alignment input
- fetch: A fetching approach to randomly subsample reads given read overlap position. Requires indexed input (.bai)
[default: stream]
--swap-distance <INT>
[Stream] A maximum distance (bp) allowed between start position of new read and the worst read in the heap to consider them to be 'swappable'.
Larger values allow swapping reads over greater distances, but may cause local undersampling. A value of `0` means only allows swap between reads that have the same start position.
[default: 5]
--step-size <INT>
[Fetch] When a region has less than the desired coverage, the step size to move along the chromosome to find more reads.
The lowest of the step and the minimum end coordinate of the reads in the region will be used. This parameter can have a significant impact on the runtime of the subsampling process.
[default: 100]
--batch-size <INT>
[Fetch] The size of the genomic window (bp) to cache into memory at once.
Larger values reduce disk seeking, but at the cost of high memory usage. The minimum value is 1,000 bp to avoid small region queries.
[default: 10000]
-@, --threads <INT>
Number of threads to use for BAM (de)compression
Only BAM benefits from this - SAM is uncompressed and CRAM has no multithreaded codec, so this is a no-op for those formats. `--threads 1` is identical to not passing this option.
[default: 1]
-h, --help
Print help (see a summary with '-h')
-V, --version
Print version
```
## Benchmark
> “Time flies like an arrow; fruit flies like a banana.”
> ― Anthony G. Oettinger
The real question is: will `rasusa` just needlessly eat away at your precious time on
earth?
To do this benchmark, I am going to use [hyperfine][hyperfine].
The data I used comes from
> [Bainomugisa, Arnold, et al. "A complete high-quality MinION nanopore assembly of an
> extensively drug-resistant Mycobacterium tuberculosis Beijing lineage strain
> identifies novel variation in repetitive PE/PPE gene regions." Microbial genomics 4.7
> (2018).][1]
> [!NOTE]
> These benchmarks are for `reads` only as there is no other tool that replicates the functionality of `aln`.
> [!NOTE]
> The tables below are regenerated automatically on each release by
> [`benches/update_readme.sh`](benches/update_readme.sh), run on a GitHub Actions
> `ubuntu-latest` runner. Absolute numbers are therefore specific to that runner and will
> vary run to run - the meaningful figure is the **relative** speedup, not the raw
> timings.
### Single long read input
Download and rename the FASTQ
```shell
URL="ftp://ftp.sra.ebi.ac.uk/vol1/fastq/SRR649/008/SRR6490088/SRR6490088_1.fastq.gz"
wget "$URL" -O - | gzip -d -c > tb.fq
```
The file size is 2.9G, and it has 379,547 reads.
We benchmark against `filtlong` using the same strategy outlined in
[Motivation](#motivation).
```shell
TB_GENOME_SIZE=4411532
COVG=50
TARGET_BASES=$(( TB_GENOME_SIZE * COVG ))
FILTLONG_CMD="filtlong --target_bases $TARGET_BASES tb.fq"
RASUSA_CMD="rasusa reads tb.fq -c $COVG -g $TB_GENOME_SIZE -s 1"
hyperfine --warmup 3 --runs 10 --export-markdown results-single.md \
"$FILTLONG_CMD" "$RASUSA_CMD"
```
#### Results
<!-- BENCH:single:START -->
| Command | Mean [s] | Min [s] | Max [s] | Relative |
|:---|---:|---:|---:|---:|
| `filtlong --target_bases 220576600 tb.fq` | 32.032 ± 0.082 | 31.922 | 32.155 | 51.95 ± 0.43 |
| `/home/runner/work/rasusa/rasusa/target/release/rasusa reads tb.fq -c 50 -g 4411532 -s 1 -o /dev/null` | 0.617 ± 0.005 | 0.608 | 0.622 | 1.00 |
**Summary**: `rasusa` ran 51.95 ± 0.43 times faster than `filtlong`.
<!-- BENCH:single:END -->
### Paired-end input
Download and then deinterleave the FASTQ with [`pyfastaq`][pyfastaq]
```shell
URL="ftp://ftp.sra.ebi.ac.uk/vol1/fastq/SRR648/008/SRR6488968/SRR6488968.fastq.gz"
wget "$URL" -O - | gzip -d -c - | fastaq deinterleave - r1.fq r2.fq
```
Each file's size is 179M and has 283,590 reads.
For this benchmark, we will use [`seqtk`][seqtk]. We will also test `seqtk`'s 2-pass
mode as this is analogous to `rasusa reads`.
```shell
NUM_READS=140000
SEQTK_CMD_1="seqtk sample -s 1 r1.fq $NUM_READS > /tmp/r1.fq; seqtk sample -s 1 r2.fq $NUM_READS > /tmp/r2.fq;"
SEQTK_CMD_2="seqtk sample -2 -s 1 r1.fq $NUM_READS > /tmp/r1.fq; seqtk sample -2 -s 1 r2.fq $NUM_READS > /tmp/r2.fq;"
RASUSA_CMD="rasusa reads r1.fq r2.fq -n $NUM_READS -s 1 -o /tmp/r1.fq -o /tmp/r2.fq"
hyperfine --warmup 10 --runs 100 --export-markdown results-paired.md \
"$SEQTK_CMD_1" "$SEQTK_CMD_2" "$RASUSA_CMD"
```
#### Results
<!-- BENCH:paired:START -->
| Command | Mean [ms] | Min [ms] | Max [ms] | Relative |
|:---|---:|---:|---:|---:|
| `seqtk sample -s 1 r1.fq 140000 > /tmp/tmp.u9aulQzBnc/o1.fq; seqtk sample -s 1 r2.fq 140000 > /tmp/tmp.u9aulQzBnc/o2.fq;` | 597.3 ± 6.9 | 578.4 | 625.6 | 1.38 ± 0.09 |
| `seqtk sample -2 -s 1 r1.fq 140000 > /tmp/tmp.u9aulQzBnc/o1.fq; seqtk sample -2 -s 1 r2.fq 140000 > /tmp/tmp.u9aulQzBnc/o2.fq;` | 456.9 ± 5.4 | 447.1 | 482.3 | 1.06 ± 0.07 |
| `/home/runner/work/rasusa/rasusa/target/release/rasusa reads r1.fq r2.fq -n 140000 -s 1 -o /tmp/tmp.u9aulQzBnc/o1.fq -o /tmp/tmp.u9aulQzBnc/o2.fq` | 432.1 ± 26.3 | 333.8 | 536.8 | 1.00 |
**Summary**: `rasusa reads` ran 1.38 times faster than `seqtk` (1-pass) and 1.06 times faster than `seqtk` (2-pass)
<!-- BENCH:paired:END -->
So, `rasusa reads` is faster than `seqtk` but doesn't require a fixed number of reads -
allowing you to avoid doing maths to determine how many reads you need to downsample to
a specific coverage. 🤓
### v4.1.0 → v5.0.0 runtime/memory improvements
<!-- INTERNAL-BENCH:START -->
> [!NOTE]
> This block is hand-written, not auto-regenerated - see [`benches/README.md`](benches/README.md)
> if you want to reproduce these numbers yourself with `benches/bench.sh compare 4.1.0 5.0.0`.
Between v4.1.0 and v5.0.0, `rasusa` went through a runtime/memory optimisation pass
(tracked in [#170](https://github.com/mbhall88/rasusa/issues/170)). The internal
benchmark harness (`benches/bench.sh`, distinct from the `filtlong`/`seqtk` comparisons
above) shows the gains per scenario, comparing `bench.sh run` on the `4.1.0` and `5.0.0`
tags:
| Scenario | Wall time | Peak RSS | Notes |
|----------------------------------------------|----------:|------------:|---------------------------------------------------|
| `reads` `--num`/`--frac` (1M reads, keep 25%) | ~-8% | ~-31/-34% | O(k) selection ([#180](https://github.com/mbhall88/rasusa/issues/180)) |
| `reads` `--num` (10M reads, keep 1,000) | ~-8% | ~-74% | same change, extreme n≫k case (54 MB → 14 MB) |
| `reads` paired `--num` | ~-7% | ~-25% | O(k) selection applies to paired mode too |
| `reads` `--coverage`/`--bases` | ~-8% | unchanged | selection algorithm here is deliberately untouched |
| `aln` (stream/fetch/paired) | ~-16/-21% | mixed | from other S1-S9 work (buffered I/O, lighter alignment domain types); RSS deltas here are within measurement noise on these tiny test fixtures |
The standout is `--num`/`--frac`'s peak memory: selection used to shuffle a full
`n`-length index vector regardless of how few reads were being kept, so memory scaled
with input size, not output size. It now scales with `k` (reads kept) instead - the
gap widens the larger the input gets relative to what you're keeping, which is exactly
the "downsample a huge run to a small target set" use case this tool exists for. See
[PR #199](https://github.com/mbhall88/rasusa/pull/199) for the full investigation and a
scaling chart across input sizes.
<!-- INTERNAL-BENCH:END -->
## Contributing
If you would like to help improve `rasusa` you are very welcome!
For changes to be accepted, they must pass the CI and coverage checks. These include:
- Code is formatted with `rustfmt`. This can be done by running `cargo fmt` in the
project directory.
- There are no compiler errors/warnings. You can check this by running `cargo clippy
--all-features --all-targets -- -D warnings`
- Code coverage has not reduced. If you want to check coverage before pushing changes, I
use [`tarpaulin`][tarpaulin].
## Citing
If you use `rasusa` in your research, it would be appreciated if you could cite it.
[](https://doi.org/10.46471/gigabyte.180)
[](https://doi.org/10.21105/joss.03941)
### Alignment subsampling (v3.0.0+)
> Furqon, A. D. C., Roberts, L. W., & Hall, M. B. (2026). Efficient downsampling of genome alignments with Rasusa. Gigabyte, 2026, gigabyte180. https://doi.org/10.46471/gigabyte.180
### General read subsampling
> Hall, M. B., (2022). Rasusa: Randomly subsample sequencing reads to a specified coverage. Journal of Open Source
> Software, 7(69), 3941, https://doi.org/10.21105/joss.03941
You can also find these citations in the [`CITATION.cff`](CITATION.cff) file or by running `rasusa cite`.
### Bibtex
You can get the following citations by running `rasusa cite`
```Bibtex
@article{furqon_efficient_2026,
title = {Efficient downsampling of genome alignments with {Rasusa}},
volume = {2026},
issn = {2709-4715},
url = {https://gigabytejournal.com/articles/180},
doi = {10.46471/gigabyte.180},
language = {en},
urldate = {2026-04-28},
journal = {Gigabyte},
publisher = {GigaScience Press},
author = {Furqon, Achmad Dimas Cahyaning and Roberts, Leah W. and Hall, Michael B.},
month = apr,
year = {2026},
pages = {gigabyte180--0}
}
@article{Hall2022,
doi = {10.21105/joss.03941},
url = {https://doi.org/10.21105/joss.03941},
year = {2022},
publisher = {The Open Journal},
volume = {7},
number = {69},
pages = {3941},
author = {Michael B. Hall},
title = {Rasusa: Randomly subsample sequencing reads to a specified coverage},
journal = {Journal of Open Source Software}
}
```
[1]: https://doi.org/10.1099/mgen.0.000188
[brew-tap]: https://github.com/brewsci/homebrew-bio
[channels]: https://bioconda.github.io/user/install.html#set-up-channels
[conda]: https://docs.conda.io/projects/conda/en/latest/user-guide/install/
[docker]: https://docs.docker.com/v17.12/install/
[dockerhub]: https://hub.docker.com/r/mbhall88/rasusa
[dpryan79]: https://github.com/dpryan79
[filtlong]: https://github.com/rrwick/Filtlong
[homebrew]: https://docs.brew.sh/Installation
[hyperfine]: https://github.com/sharkdp/hyperfine
[tarpaulin]: https://github.com/xd009642/tarpaulin
[log-lvl]: https://docs.rs/log/0.4.6/log/enum.Level.html#variants
[mgen-ref]: https://doi.org/10.1099/mgen.0.000294
[one-pass-accuracy]: benches/README.md#one-pass-fraction-sampling-how-close-does-it-land
[pr-help]: https://github.com/bioconda/bioconda-recipes/pull/18690
[pyfastaq]: https://github.com/sanger-pathogens/Fastaq
[quay.io]: https://quay.io/repository/mbhall88/rasusa
[rust]: https://www.rust-lang.org/tools/install
[score]: https://github.com/rrwick/Filtlong#read-scoring
[seed]: https://en.wikipedia.org/wiki/Random_seed
[seqtk]: https://github.com/lh3/seqtk
[apptainer]: https://github.com/apptainer/apptainer
[snakemake]: https://snakemake.readthedocs.io/en/stable/
[triples]: https://clang.llvm.org/docs/CrossCompilation.html#target-triple
[wrapper]: https://snakemake-wrappers.readthedocs.io/en/stable/wrappers/rasusa.html
[zstd]: https://github.com/facebook/zstd
[ghcr]: https://github.com/mbhall88/rasusa/pkgs/container/rasusa