Phlya/tagmap
Map insertions of transopable elements using tagmentation
Overview
Latest release: None, Last update: 2026-08-12
Share link: https://snakemake.github.io/snakemake-workflow-catalog?wf=Phlya/tagmap
Quality control: linting: passed formatting: failed
Wrappers: bio/bwa-memx/index bio/bwa-memx/mem bio/fastp bio/fastqc
Workflow Rule Graph
This visualization of the workflow’s rule graph was automatically generated using Snakevision
Deployment
Step 1: Install Snakemake and Snakedeploy
Snakemake and Snakedeploy are best installed via the Conda package manager. It is recommended to install conda via Miniforge. Run
conda create -c conda-forge -c bioconda -c nodefaults --name snakemake snakemake snakedeploy
to install both Snakemake and Snakedeploy in an isolated environment. For all following commands ensure that this environment is activated via
conda activate snakemake
For other installation methods, refer to the Snakemake and Snakedeploy documentation.
Step 2: Deploy workflow
With Snakemake and Snakedeploy installed, the workflow can be deployed as follows. First, create an appropriate project working directory on your system and enter it:
mkdir -p path/to/project-workdir
cd path/to/project-workdir
In all following steps, we will assume that you are inside of that directory. Then run
snakedeploy deploy-workflow https://github.com/Phlya/tagmap . --tag None
Snakedeploy will create two folders, workflow and config. The former contains the deployment of the chosen workflow as a Snakemake module, the latter contains configuration files which will be modified in the next step in order to configure the workflow to your needs.
Step 3: Configure workflow
To configure the workflow, adapt config/config.yml to your needs following the instructions below.
Step 4: Run workflow
The deployment method is controlled using the --software-deployment-method (short --sdm) argument.
To run the workflow with automatic deployment of all required software via conda/mamba, use
snakemake --cores all --sdm conda
Snakemake will automatically detect the main Snakefile in the workflow subfolder and execute the workflow module that has been defined by the deployment in step 2.
For further options such as cluster and cloud execution, see the docs.
Step 5: Generate report
After finalizing your data analysis, you can automatically generate an interactive visual HTML report for inspection of results together with parameters and code inside of the browser using
snakemake --report report.zip
Configuration
The following section is imported from the workflow’s config/README.md.
Configuration
Copy example_config.yaml to config.yaml and edit it for your experiment.
Every key is documented inline; only refgen_path, chrom_sizes_path,
cassette_name, forward_primer_sequence and reverse_primer_sequence are
required, everything else has a default (see
workflow/schemas/config.schema.yaml for the authoritative list, which
snakemake validates the config against on every run).
Sample sheets
Give either or both, depending on what you’re analysing:
samples.tsv(samples_path) - paired-end NGS tagmentation libraries. Columns:name,fastq1,fastq2. Rows sharing anameare merged, so a library split across lanes can be listed as several rows.sanger_samples.tsv(sanger_samples_path) - Sanger runs of individual long reads, for pinpointing the integration site in a clonal line. Columns:name, and eitherab1(a folder or glob of.ab1traces) orfastq(already-basecalled reads); optionallydirection(forward/reverse),clone, andngs_sample(an NGS sample fromsamples.tsvcovering the same material, used to cross-validate the Sanger site - seevalidate_sanger_with_ngs).
Read-level metadata such as well, clone and primer direction is normally
pulled from the trace file names via sanger_name_regex; sample-sheet columns
fill in whatever the regex doesn’t capture, or override it for a whole row.
Both sheets take #-prefixed comment lines and blank rows.
Output paths
Every stage’s output folder (fastq_folder, bams_folder, pairs_folder, …)
and the primer_position_file/fasta_index_file can be set individually, but
default to living under results_folder (itself defaulting to results), so
usually just setting results_folder is enough to get a self-contained,
namespaced results tree - only override an individual path if you want that
one file somewhere else.
Summary tables
Besides the site calls, the workflow writes small TSV summaries to
stats_folder, and assembles them into a single report.md:
ngs_qc_stats.tsv- per NGS library: read pairs anchored at an ITR primer (evidence of mobilization), and how many of the resulting insertion sites were seen from both sides of the cassette versus only one.sanger_qc_stats.tsv- Sanger reads per run and primer direction, how many passed QC, and a summary of why the rest didn’t.sanger_clone_summary.tsv- per clone, whether the forward and reverse primers each gave a confirmed site, failed QC, or were never sequenced, andwhen NGS data for the same material is available - whether that clone’s site was independently confirmed there, and from which side(s). When
original_insertion_site/original_insertion_upstream_seqis configured, a side whose site lands there instead of a new locus is classifiedunmobilizedrather thanconfirmed, and the clone gets anunmobilizedcolumn.
sanger_positions.tsv- distinct genomic positions found across the clustered Sanger sites, per plate and in total (a locus seen on more than one plate, e.g. an unmobilized founder control resequenced, still counts once).validation_summary.tsv- Sanger sites confirmed by the NGS data, per run and in total (only written whenvalidate_sanger_with_ngsapplies).
Workflow parameters
The following table is automatically parsed from the workflow’s config.schema.y(a)ml file.
Parameter |
Type |
Description |
Required |
Default |
|---|---|---|---|---|
samples_path |
string |
TSV of paired-end NGS libraries (columns name, fastq1, fastq2). Rows sharing a name are merged. Leave unset for a Sanger-only run. |
| | |
| sanger_samples_path | string | TSV of Sanger runs (columns name, ab1 or fastq, and optionally direction, clone, ngs_sample). Leave unset for an NGS-only run.
| | |
| refgen_path | string | Reference FASTA, including the cassette as an extra contig. | yes | |
| chrom_sizes_path | string | Chrom sizes of refgen_path, including the cassette contig. | yes | |
| chrom_sizes_path_no_cassette | string | Chrom sizes of the genome alone, without the cassette. Used to trim output for genome browsers, which do not know about the cassette.
| | |
| mapper | string | From slowest/leanest to fastest/most memory hungry. minibwa (https://github.com/lh3/minibwa) is a newer, faster aligner that isn’t supported by the bwa-memx wrapper the other three use, so it’s mapped with its own rules: minibwa mem for the single-end Sanger reads, minibwa map –hic for paired-end NGS reads (matching the -SP flags used for the other mappers - independent mate mapping, no proper-pair assumption, since these are chimeric tagmentation fragments rather than conventional inserts).
| | bwa-mem2 |
| cassette_name | string | Name of the cassette contig in refgen_path. | yes | |
| forward_primer_sequence | string | The part of the forward primer that binds the ITR. If several PCRs were used, the primer of the last one. Used for the NGS branch, and as the Sanger branch’s default unless sanger_forward_primer_sequence is set. Required when samples_path is set; a Sanger-only project may instead give just sanger_forward_primer_sequence/sanger_reverse_primer_sequence.
| | |
| reverse_primer_sequence | string | As forward_primer_sequence, for the other end of the cassette. | | |
| sanger_forward_primer_sequence | string | Overrides forward_primer_sequence for the Sanger branch only, for projects where the actual Sanger sequencing primer differs from the primer used to build the NGS library (e.g. a different primer against the same ITR), or where a Sanger-only project has no NGS primer to default to at all. Defaults to forward_primer_sequence; required (along with sanger_reverse_primer_sequence) when sanger_samples_path is set and forward_primer_sequence/reverse_primer_sequence are not given.
| | |
| sanger_reverse_primer_sequence | string | As sanger_forward_primer_sequence, for the reverse primer. | | |
| insertion_seq | string | Motif the transposon integrates into - TA for Sleeping Beauty, TTAA for PiggyBac. Sites are snapped onto the nearest occurrence.
| | TA |
| do_fastqc | boolean | | | false |
| trim | boolean | Trim reads with fastp before mapping. | | false |
| trim_args | string | Extra arguments passed to fastp when trim is true. | | |
| dedup | boolean | | | true |
| blacklist | string | Peaks overlapping this BED are dropped. Unset to disable. | | |
| peak_caller | string | Which pairs are used, and how, in one setting.
A TagMap fragment runs [tagmentation site] –genome– [insertion site] –cassette– [ITR primer]. When a read is long enough to carry on past the junction, the first genomic base it reaches is the insertion site itself, the same base for every molecule. When it is not, the only genomic coordinate in the pair comes from the other read, which starts at the tagmentation cut - a different random position per molecule.
“coverage” clusters those genomic coordinates and cuts on read count, without distinguishing the two. Where reads rarely span the junction that makes one peak per molecule: on the 2G6 MiSeq run, 30,053 sites for a pool of ~800.
“coverage_junctions” first throws away every pair that did not read through, then does the same clustering with an automatic intensity threshold. Positions are exact, but everything the shorter reads saw is gone, and the threshold drops rare sites.
“junction_tiered” keeps both kinds and uses each for what it is: junctions place a site to the base and are believed on their own, ITR-anchored pairs may only corroborate one, or form a peak that must still clear min_peak_reads/min_peak_frac. It also counts distinct tagmentation positions rather than reads, so abundance is captured molecules instead of PCR. It is the default, and needs no up-front decision about read length.
| | junction_tiered |
| junction_jitter | integer | junction_tiered only. Junction positions this close apart are one site, absorbing mapping wobble. Keep it small: in a local re-mobilisation assay, distinct real insertions can sit only 1-2bp apart and still show cleanly opposite orientations read-for-read - set to 0 to trust every distinct junction base as its own site (still subject to max_dist_between_sides for combining a site’s forward/reverse peaks).
| | 2 |
| max_fragment | integer | junction_tiered only. Longest junction-to-tagmentation distance a pair may span, i.e. the largest library fragment. Bounds how far an ITR-anchored pair may reach to corroborate a junction.
| | 1000 |
| min_junction_frags | integer | junction_tiered only. Distinct tagmentation positions among junction-spanning reads needed to call a site. The default of 1 is the point of the tiering - one read that spelled out the junction is proof of an integration, however rare it is in the pool.
| | 1 |
| min_junction_mapq | integer | junction_tiered only. Junction pairs below this MAPQ are dropped rather than trusted on their own. Worth raising above 0 if single-read sites in repeats are a concern.
| | 0 |
| motif_max_support | integer | junction_tiered only. Junction peaks supported by this many distinct molecules or fewer must also carry insertion_seq at the junction, which is what a genuine integration looks like. 0 disables the check. It is what makes min_junction_frags: 1 safe: on the 2G6 long-read library a singleton junction carries a TA only 48% of the time against a 23% background - so barely a third can be real - while junctions seen on four or more molecules reach 90%. Asking thin peaks to look like integrations keeps the genuinely rare ones without the chimeras, and calibrates itself, since the motif rate against background is a false-discovery estimate needing no labelled sites.
| | 0 |
| min_peak_width | integer | | | 1 |
| min_peak_reads | integer | | | 10 |
| min_peak_frac | number | | | 0.01 |
| min_peak_dist | integer | Same-side peak positions within this distance are merged into one peak. Suited to a genome-wide screen, where independent insertions land far apart and anything closer is read scatter around one true site. Turn this down (e.g. to ~5) for an assay that packs many independent insertions into one small locus (e.g. local re-mobilization of a cassette), where true neighbouring sites can be under 100bp apart - left at the default, such sites get silently merged and undercounted.
| | 1000 |
| min_peak_positions | integer | | | 1 |
| snap_window | integer | Before forward/reverse peaks are combined, each one is snapped to the nearest insertion_seq motif within this many bases, so both sides of one real integration land on the exact same base regardless of a few bp of slop (the +/- strand pos5 reporting offset, indel wobble, Sleeping Beauty’s own target-site duplication) - checked against real data, a peak’s own raw position can be off by up to ~4bp even when correct. This is what lets max_dist_between_sides below be set tight: with snapping on, two peaks from one real integration should already share a position, so a small max_dist is only covering the cases snapping itself could not resolve (no motif within window). 0 disables snapping.
| | 5 |
| max_dist_between_sides | integer | Peaks from the forward and reverse side of the cassette within this distance are combined into one insertion site. With snap_window enabled (the default), this only needs to be small - e.g. 2 - since real pairs should already coincide; it still governs the genome-wide- vs-dense-locus tradeoff (see min_peak_dist) for whatever snapping could not resolve.
| | 100 |
| min_orientation_support | integer | A peak needs at least this many distinct molecules before it gets a say in which way its site faces. Checked against Sanger-confirmed sites: peaks under ~20 molecules were barely better than a coin flip on strand, while every peak above ~100 molecules was correct. Raise it to demand more evidence before committing to a strand (more sites fall back to “.”); 0 disables the check and trusts every peak’s own orientation, however thin.
| | 20 |
| sanger_basecaller | string | abi-trim reads the base calls embedded in the .ab1 trace and applies Mott trimming; tracy re-calls bases from the trace itself.
| | abi-trim |
| sanger_map_args | string | Extra arguments for mapping Sanger reads. | | -T 30 -k 19 |
| sanger_min_mapq | integer | Minimum MAPQ for a genomic alignment to define a site. | | 30 |
| sanger_min_aligned | integer | Minimum aligned length of the genomic segment of a read. | | 30 |
| sanger_max_unexplained | integer | A read fails QC if more than this many bases before the genomic segment are neither aligned to the genome nor to the cassette.
| | 50 |
| sanger_name_regex | [‘string’, ‘null’] | Regex with named groups matched against each .ab1/read file stem, e.g. ‘(?P
Linting and formatting
Linting results
All tests passed!
Formatting results
1[DEBUG]
2[DEBUG] In file "/tmp/tmpf_g9idw3/workflow/rules/report.smk": Formatted content is different from original
3[DEBUG]
4[DEBUG]
5[DEBUG]
6[DEBUG]
7[DEBUG]
8[DEBUG] In file "/tmp/tmpf_g9idw3/workflow/rules/common.smk": Formatted content is different from original
9[DEBUG]
10[DEBUG]
11[INFO] 2 file(s) would be changed 😬
12[INFO] 6 file(s) would be left unchanged 🎉
13
14snakefmt version: 0.11.5