sourmash plugin for repeat-robust mutation rate estimation (r_pp, r_pc, r_cc).
Project description
Repeat-robust estimators for mutation rate
This repo describes a sourmash plugin that estimates the mutation rate between two DNA sequences. It implements the estimators from this ISMB paper. These estimators are well-suited for highly repetitive sequences such as centromeres.
Installation
# Option 1: conda + pip (recommended)
conda install -c conda-forge -c bioconda sourmash
pip install sourmash-plugin-repeat-robust-mutation-rate-estimators
# Option 2: pip only
pip install sourmash
pip install sourmash-plugin-repeat-robust-mutation-rate-estimators
Verify the plugin is recognized:
sourmash scripts
You should see sketch and mutation_rate listed under available plugin commands.
Quick start
To start, you need two FASTA files (.fa or .fasta) as input. Each FASTA file may contain multiple records. This tool is designed for full-length assembled sequences. We currently do not support sequencing reads as input.
You first need to sketch each of the FASTA files. There are three sketch modes: standard, multipliticy, and extended. We will describe which one to use later. As an example, if your two files are s.fa and t.fa and you want to sketch them using standard and multiplicity modes respectively, then you can run:
sourmash scripts sketch s.fa --sketch-mode standard -o s.sig -k 21 --scaled 1000
sourmash scripts sketch t.fa --sketch-mode multiplicity -o t.sig -k 21 --scaled 1000
This uses a k-mer size of 21 and a subsampling rate of 1/1000. A smaller subsampling rate leads to a smaller sketch size and faster downstream estimation, but introduces more variance in downstream estimation.
Once you've sketched the files, you are ready to compare them. There various estimators which we describe below. For example, to compare them using the r_pc estimator use
sourmash scripts mutation_rate --estimator r_pc --s-sig s.sig --t-sig t.sig
This will generate output that looks like
Estimator : r_pc
k : 21
scaled : 1000
L_s : 4800000
Estimated mutation rate : 0.012345
Choice of $s$ and $t$
Our estimators are not symmetric and your decision on which fasta file to label as $s$ and which as $t$ may give different results. Our underlying model assumes that there is a starting string $s$ where each position has a substitution at rate $r$, resulting in the mutated string $t$. To choose which of your FASTA files to use as $s$ and which as $t$, we give the following guidelines:
- If you know the biological direction (e.g., one sequence is ancestral), assign the ancestral sequence as $s$ and the derived sequence as $t$.
- If you are unsure, use the longer sequence as $s$ and the shorter as $t$.
- If you are comparing many sequences (e.g., building a phylogenetic tree), be consistent in the way you assign $s$ and $t$, e.g. always assign the longer sequence as $s$ across all pairs.
- If the biological direction of the sequences is unclear, we recommend running the estimator twice, swapping the roles of $s$ and $t$ between the runs. That way, you can compare the difference in the estimates to the uncertainty you are willing to tolerate.
Choice of estimator
We provide three estimators, and each estimator requires a specific sketch mode for $s$ and $t$.
| Estimator | Sketch mode for s | Sketch mode for t | Accuracy |
|---|---|---|---|
r_pp |
standard |
standard |
Lower |
r_pc |
standard |
multiplicity |
Medium |
r_cc |
extended |
multiplicity |
Highest |
r_pp(presence-to-presence): ignores the multiplicity of $k$-mers of both $s$ and $t$. Use this estimator if you do not trust the multiplicity information in your input sequences, e.g. you expect collapsed repeats or the input is unitigs.r_pc(presence-to-count): ignores the multiplicity of $k$-mers in $s$ but uses the multiplicity of $k$-mers in $t$.r_cc(count-to-count): uses $k$-mer multiplicites of both $s$ and $t$, with a bias correction term precomputed from $s$. This is the most accurate estimator, especially for highly repetitive sequences. However, sketching $s$ inextendedmode will take longer.
To estimate the mutation rate, you can use one of the following commands:
sourmash scripts mutation_rate --estimator r_pp --s-sig s.sig --t-sig t.sig
sourmash scripts mutation_rate --estimator r_pc --s-sig s.sig --t-sig t.sig
sourmash scripts mutation_rate --estimator r_cc --s-sig s.sig --t-sig t.sig
If the estimator finds that the sketch modes do not match its requirements, it will report an error and tell you which modes are needed.
Support
Please file issues at https://github.com/Wu-Haonan/sourmash-plugin-repeat-robust-mutation-rate-estimators/issues
Citation
If you use this plugin, please cite:
Wu, H. and Medvedev, P. (2026). The gift of novelty: repeat-robust k-mer-based estimators of mutation rates. Proceedings of ISMB 2026, to appear in Bionformatics.
License
MIT License. See LICENSE for details.
Project details
Release history Release notifications | RSS feed
Download files
Download the file for your platform. If you're not sure which to choose, learn more about installing packages.
Source Distribution
Built Distribution
Filter files by name, interpreter, ABI, and platform.
If you're not sure about the file name format, learn more about wheel file names.
Copy a direct link to the current filters
File details
Details for the file sourmash_plugin_repeat_robust_mutation_rate_estimators-0.1.2.tar.gz.
File metadata
- Download URL: sourmash_plugin_repeat_robust_mutation_rate_estimators-0.1.2.tar.gz
- Upload date:
- Size: 9.1 kB
- Tags: Source
- Uploaded using Trusted Publishing? No
- Uploaded via: twine/6.2.0 CPython/3.12.13
File hashes
| Algorithm | Hash digest | |
|---|---|---|
| SHA256 |
6ed603a9f85dae325b67293a5f698e847b60f7129f4bae0fbc5860566a8c7db2
|
|
| MD5 |
d57c29209cd85aed33e197dd622e29da
|
|
| BLAKE2b-256 |
aec4f35a875448daafdcae0e9be7b36bfcdce18a79621126cdab08b8b7e7cb23
|
File details
Details for the file sourmash_plugin_repeat_robust_mutation_rate_estimators-0.1.2-py3-none-any.whl.
File metadata
- Download URL: sourmash_plugin_repeat_robust_mutation_rate_estimators-0.1.2-py3-none-any.whl
- Upload date:
- Size: 8.7 kB
- Tags: Python 3
- Uploaded using Trusted Publishing? No
- Uploaded via: twine/6.2.0 CPython/3.12.13
File hashes
| Algorithm | Hash digest | |
|---|---|---|
| SHA256 |
38573181b25aac95bf647823f4bcaeafe9a8167e6f2b9c4f6a6e48b228b6b448
|
|
| MD5 |
012c96289f53415508270874e5272d52
|
|
| BLAKE2b-256 |
931ee4b7aa222ecf38a76f4e0c12970bcc800027a60bed903c894e69eaa6b53a
|