This Nextflow DSL2 workflow demultiplexes paired-end bulk nano-CUT&Tag reads
with an 8-base I2 barcode, aligns each derived sample, calls matched-IgG
NanoScope-compatible broad peaks, calculates library/peak QC, and optionally
tests expected and de novo motifs or peak-overlap enrichment against public
ChIP-seq references. The repository includes the exact six-library manifest
mapping in assets/samples.example.csv.
- A POSIX shell and Java suitable for the installed Nextflow release.
- Nextflow 23.10.0 or later (enforced by
manifest.nextflowVersion). - One supported dependency backend:
- Conda, Mamba, or Micromamba with
-profile conda; or - Docker with
-profile docker.
- Conda, Mamba, or Micromamba with
- A reference FASTA or complete Bowtie2 index.
- The input manifest and readable R1, R2, and I2 FASTQs.
- A MACS2 genome-size shortcut such as
hs, or a positive effective genome size integer. - Python 3.10 or newer on the launch host when
--chipseq_inputis used. Before workflow scheduling, this preflight streams FASTA chromosome lengths and runs the repository CSV/BED validator. Enrichment tasks themselves still use the selected Conda or Docker profile.
The conda profile is the most direct portable installation path. Nextflow
automatically creates each process environment from the repository-local,
exactly pinned YAML file when it is first needed; users should not create a
single combined environment by hand. The Docker profile uses the
process-specific pinned images declared in the modules.
Nextflow runtime execution was not verified in this workspace. Unit tests, static assertions, and direct synthetic fixtures were verified, but Nextflow was unavailable; therefore DSL2 execution, Conda solving, container pulls, real bioinformatics-tool compatibility, and cache/resume behavior still need verification on a supported runtime.
The canonical environments are:
| File | Direct pinned tool | Used for |
|---|---|---|
envs/python.yml |
Python 3.12.3, matplotlib 3.9.2 | Manifest validation, demultiplexing, enrichment statistics/plots, custom QC, pipeline metadata |
envs/fastqc.yml |
FastQC 0.12.1 | Read QC |
envs/bowtie2.yml |
Bowtie2 2.5.4 | Index construction and alignment |
envs/samtools.yml |
SAMtools 1.20 | BAM sorting/filtering/indexing, metrics, fragments |
envs/deeptools.yml |
deepTools 3.5.5 | RPKM bigWig and optional TSS enrichment |
envs/macs2.yml |
MACS2 2.2.9.1 | Broad primary and narrow motif-only peak calls |
envs/bedtools.yml |
BEDTools 2.31.1 | Blacklist filtering and motif sequence/background intervals |
envs/meme.yml |
MEME Suite 5.5.7 | AME, STREME, and FIMO |
envs/multiqc.yml |
MultiQC 1.25.2 | Consolidated report and summary tables |
Each environment uses conda-forge before bioconda. On a local workstation,
the automatic cache defaults to .nextflow-conda under the launch directory:
nextflow run main.nf \
-profile conda \
--input assets/samples.example.csv \
--fasta /references/GRCh38.fa \
--macs_genome_size hs \
--outdir resultsFor a cluster, place both the Conda cache and work directory on storage visible with the same path from every compute node:
export NXF_CONDA_CACHEDIR=/shared/nanocut/nextflow-conda
export NXF_HOME=/shared/nanocut/nextflow-home
nextflow run main.nf \
-profile conda \
-work-dir /shared/nanocut/work \
--input /shared/project/samples.csv \
--fasta /shared/references/GRCh38.fa \
--macs_genome_size hs \
--outdir /shared/project/resultsNXF_CONDA_CACHEDIR overrides conda.cacheDir; without it the local automatic
cache is used. NXF_HOME is optional and controls Nextflow's own home/cache
location, not the pipeline result directory. The shared directories must be
writable during initial environment creation and readable by all jobs.
Mamba/Micromamba selection may be supplied in a user Nextflow configuration;
the repository YAML files remain the dependency source of truth.
The offline test profile modifies PATH, LC_ALL, LANG, and
PYTHONDONTWRITEBYTECODE only inside its processes so repository fake tools
are deterministic. Those are test-only environment variables and are not
production configuration knobs.
The CSV has one row per barcode-derived sample and this exact header:
sample_id,library_id,input_group,barcode,assay_target,is_control,control_id,expected_motif,r1,r2,i2The columns mean:
sample_id: unique derived-sample identifier. It must begin with a letter or digit and contain only letters, digits,.,_, or-.library_id: physical-library identifier, with the same character rule. Rows from one library repeat the same FASTQ paths intentionally; the physical library is streamed only once.input_group: biological input group used to require same-input matched controls, for example25K,50K, or100K.barcode: expected I2 sequence in FASTQ orientation. All manifest barcodes must be non-empty and equal length, and barcodes within a library must be unique.assay_target: one ofIgG,CTCF,GATA1, orRUNX1.is_control:true,1, oryesfor controls;false,0, ornofor targets, case-insensitively.control_id: matched IgGsample_id. It is blank for controls and required for targets. The referenced row must be an IgG control with the sameinput_group; it may come from another physical library.expected_motif: regular expression matched against motif ID and alternate name. It is blank for controls and required for every target.r1,r2,i2: synchronized FASTQ paths. Plain or gzip-compressed content is accepted. Relative paths resolve from the manifest directory, and all three paths must be identical across rows sharing alibrary_id.
Validation stops before expensive work for missing/blank required fields,
unsupported targets, invalid identifiers/booleans, missing FASTQs, duplicate
sample IDs, duplicate within-library barcodes, unequal barcode lengths,
inconsistent library FASTQ paths, invalid controls, or cross-input control
links. R1/R2/I2 record counts and normalized read identifiers are also checked
during streaming demultiplexing. A unique closest barcode within
--barcode_mismatches is assigned; ties are ambiguous and out-of-threshold
reads are unassigned.
| Physical library | Input | TATAGCCT (barcode A) |
ATAGAGGC (barcode B) |
Matched IgG for targets |
|---|---|---|---|---|
| NX702 | 25K | NX702_IgG | NX702_CTCF | NX702_IgG |
| NX703 | 50K | NX703_IgG | NX703_CTCF | NX703_IgG |
| NX701 | 100K | NX701_IgG | NX701_CTCF | NX701_IgG |
| NX704 | 25K | NX704_GATA1 | NX704_RUNX1 | NX702_IgG |
| NX705 | 50K | NX705_GATA1 | NX705_RUNX1 | NX703_IgG |
| NX706 | 100K | NX706_GATA1 | NX706_RUNX1 | NX701_IgG |
Thus NX701–NX703 use barcode A for IgG and barcode B for CTCF, while
NX704–NX706 use barcode A for GATA1 and barcode B for RUNX1. Copy
assets/samples.example.csv, update its FASTQ paths, and retain the listed
control relationships. The expected motif expressions are CTCF, GATA1,
and RUNX1 for their respective targets.
Use one of:
--fasta /path/genome.fa: a readable FASTA whose first nonblank line is a>header. The pipeline builds and publishes a Bowtie2 index.--bowtie2_index /path/index/prefix: a prefix for exactly one complete six-file.bt2or.bt2lset:1,2,3,4,rev.1, andrev.2. Supply the prefix only, without a suffix. A supplied index takes precedence for alignment.
--fasta is still required when --motif_db or --chipseq_input is supplied
because motif windows and enrichment backgrounds use the reference sequence.
Optional --blacklist is a BED file removed from final broad peaks, motif
inputs, and enrichment backgrounds. Optional --tss_bed must be strand-aware
BED6 and takes precedence over --gtf; otherwise transcript features in the
GTF are converted to strand-aware single-base TSS records.
Conda:
nextflow run main.nf \
-profile conda \
--input samples.csv \
--fasta GRCh38.fa \
--outdir results \
--blacklist hg38-blacklist.v2.bed \
--gtf gencode.annotation.gtf \
--motif_db JASPAR2026_CORE_vertebrates.meme \
--macs_genome_size hsDocker:
nextflow run main.nf \
-profile docker \
--input samples.csv \
--bowtie2_index /references/grch38/bowtie2/grch38 \
--outdir results \
--macs_genome_size 2913022398The second example omits motif analysis; add both --fasta and --motif_db
to enable it. A user config can override executor/resources without modifying
the pipeline:
nextflow run main.nf \
-profile conda \
-c cluster.config \
--input samples.csv \
--fasta GRCh38.fa \
--macs_genome_size hs \
--outdir resultsThe deterministic offline orchestration fixture is:
nextflow run main.nf -profile testIt uses tiny bundled data and fake command implementations; it is not a biological validation and was not run here because Nextflow was unavailable.
Every pipeline CLI parameter is listed below. Nextflow runtime options such as
-profile, -c, -work-dir, and -resume use one leading dash; pipeline
parameters use two.
| Parameter | Required/default | Meaning |
|---|---|---|
--input |
Required | CSV manifest, one row per barcode-derived sample. |
--fasta |
Required unless --bowtie2_index is supplied |
Reference FASTA; also required for motif analysis. |
--bowtie2_index |
Required unless --fasta is supplied |
Complete Bowtie2 index prefix; takes precedence for alignment. |
--outdir |
results |
Result directory. |
--blacklist |
Absent | Optional BED regions removed from final broad peaks and motif inputs. |
--gtf |
Absent | Optional GTF for transcript-derived, strand-aware TSS positions. |
--tss_bed |
Absent | Optional BED6 TSS file; takes precedence over --gtf. |
--motif_db |
Absent | MEME-format database beginning with MEME version and containing at least one MOTIF record; enables motif analysis. |
--chipseq_input |
Absent | Public ChIP-seq CSV with reference_id,tf,peak_file; enables peak-overlap enrichment and requires --fasta. |
--barcode_mismatches |
0 |
Non-negative I2 Hamming-distance threshold, smaller than barcode length. |
--allow_empty |
false |
Permit a derived sample with zero assigned read pairs for diagnostic runs. |
--min_mapq |
5 |
Integer 0–255 used for the filtered BAM, coverage, and fragment QC. |
--macs_genome_size |
Required | Positive effective genome size or path-safe MACS2 shortcut such as hs. |
--macs_llocal |
Fixed 100000 |
NanoScope-compatible MACS2 local lambda window; other values are rejected. |
--macs_keep_dup |
Fixed 1 |
NanoScope-compatible MACS2 duplicate setting; other values are rejected. |
--macs_broad_cutoff |
Fixed 0.1 |
Primary broad-peak cutoff; other values are rejected. |
--macs_max_gap |
Fixed 1000 |
Primary broad-peak maximum gap; other values are rejected. |
--motif_use_narrow_peaks |
true |
Use matched-control narrow summits for motif windows. If false, use final broad-peak midpoints. |
--motif_window |
200 |
Positive total motif-window width in bp. Windows shift at reference edges; a shorter contig yields its full length. |
--enrichment_permutations |
1000 |
Positive number of attempted null permutations for every foreground/reference/model comparison. |
--enrichment_seed |
1729 |
Integer base seed recorded in enrichment outputs; comparison/model-specific streams remain deterministic. |
--enrichment_gc_tolerance |
0.02 |
GC fallback tolerance in (0, 1]; the sampler first tries a 1-percentage-point match, then this tolerance. |
Boolean CLI values must reach Nextflow as booleans (true or false).
Optional features are skipped only when their path parameter is absent; a
supplied invalid path fails validation.
The principal layout under --outdir is:
results/
demultiplex/<library_id>/
fastqc/<sample_id>/
alignment/<sample_id>/
coverage/<sample_id>/
peaks/<sample_id>/broad/raw/
peaks/<sample_id>/broad/final/
peaks/<sample_id>/narrow_motif_qc/
qc/library/<sample_id>/
qc/fragments/<sample_id>/
qc/peaks/<sample_id>/
qc/tss/<sample_id>/
motifs/<sample_id>/
enrichment/
peak_enrichment.tsv
enrichment_status.tsv
plots/
reports/multiqc/
reports/qc_dashboard/
reports/summary/
pipeline_info/
demultiplex/contains per-sample paired FASTQs plus library JSON/TSV assignment metrics.fastqc/contains paired FastQC HTML/ZIP outputs.alignment/contains the primary analysis BAM/BAI, MAPQ-filtered BAM/BAI, Bowtie2 summary, and process version records.coverage/contains RPKMcoverage.RPKM.bwfiles generated with MAPQ filtering, 50-bp bins, centered/extended reads, 250-bp smoothing, and duplicate ignoring.peaks/retains the matched-IgG raw broadPeak/gappedPeak/XLS/log, the blacklist-filteredfinal.broadPeak, and optional motif-only narrowPeak and summits.qc/library/contains SAMtools flagstat/stats/idxstats, insert sizes, markdup metrics, filtered-read/fragment counts, and custom summary inputs.qc/fragments/andqc/peaks/contain target-only BEDPE fragments, fragment-based FRiP, peak metrics, width histograms, and per-peak counts.qc/tss/contains optional TSS BED, matrix, table, profile, and status.motifs/contains foreground/background intervals and FASTAs, AME known enrichment, STREME de novo discovery, FIMO scans, status/log files, andexpected_motif_qcJSON/TSV/position outputs.enrichment/contains the complete long-formpeak_enrichment.tsv, status table, four TSV matrices and heatmap PNGs, andobserved_vs_null.png.reports/multiqc/multiqc_report.htmlis the top-level report.reports/summary/combined_target_qc.tsvis the combined target broad-peak summary.reports/multiqc/multiqc_report.htmlremains the MultiQC report.reports/qc_dashboard/qc_dashboard.htmlis the detailed, self-contained consolidated QC dashboard. It has no externalhttps://dependencies and is safe to copy with a results directory.reports/qc_dashboard/qc_summary.tsvandqc_summary.jsonare reusable, joined per-sample data products.top_motifs.tsvis the reusable top-ten AME motif summary andtss_profiles.tsvis the reusable tidy TSS profile data product. The TSS scalar is inqc_summary.tsvandqc_summary.json.reports/summary/combined_target_qc.tsvremains the combined target broad-peak summary.pipeline_info/contains the normalized manifest, validated parameter JSON, software versions, completion summary, built index when applicable, execution report, timeline, trace, and DAG.
Empty target peak sets are valid outputs with zero/NA QC and recorded motif skip status. IgG libraries receive read/alignment/library QC but are not peak-called against themselves and do not receive target FRiP by default.
Supply --chipseq_input to compare every final non-control CUT&Tag target
against public ChIP-seq peaks. The CSV contract is:
reference_id,tf,peak_file
ctcf_encode,CTCF,references/ctcf.bedreference_id must be unique and path-safe, tf is a display label, and
peak_file is a readable BED-like file with zero-based, half-open intervals.
Relative peak paths resolve from the CSV directory. Quoted CSV values are
supported. Extra public columns are ignored: public rows are always normalized
as reference_type=chipseq; reference_type is reserved for the pipeline's
internal manifest, where same-run target sets are tagged called_tf. Invalid
identifiers, duplicate IDs, malformed intervals, unknown chromosomes, and
missing files fail before enrichment is scheduled.
When at least two called TFs are present, each is also used as a reference for the others; self-comparisons are excluded. Each foreground/reference pair gets four independent null models:
random: samples chromosome and width from the foreground distributions.length_matched: preserves each foreground interval's chromosome and exact width.gc_matched: preserves each interval's chromosome, samples widths from the foreground distribution, and matches GC content.length_gc_matched: preserves chromosome and exact width and matches GC.
Both GC-aware models try a 1-percentage-point match before falling back to
--enrichment_gc_tolerance (2 percentage points by default). Every model
excludes blacklist intervals and within-replicate overlaps.
The long-form result schema is:
foreground_id,foreground_tf,reference_id,reference_type,reference_tf,
background_model,foreground_peak_count,reference_peak_count,
observed_overlap_count,null_mean_overlap,null_sd_overlap,enrichment_ratio,
empirical_p_value,permutations_requested,permutations_succeeded,seed,status
status=ok requires every requested permutation to succeed and a nonzero null
mean. insufficient_background means one or more permutations failed and has
blank inferential statistics; zero_null_mean explicitly marks a completed
null distribution whose ratio and p-value are not reported.
no_foreground_peaks and no_reference_peaks describe empty inputs. Do not
interpret any non-ok row as evidence of enrichment.
Canonical artifacts are published under enrichment/. MultiQC adds a table
only for ok rows, keys each row by the foreground/reference/model
comparison, and embeds all four configured heatmap images plus the
observed-versus-null plot. The completion summary links the canonical
peak_enrichment.tsv, observed_vs_null.png, and each heatmap when enrichment
is present; no enrichment section or links are created when the option is
absent.
These are overlap statistics, not proof of direct binding or regulatory causality. Reference cell type, assay quality, genome build, peak-calling choices, blacklist coverage, and peak-set size can dominate the result. Empirical p-value resolution is limited by the permutation count, and the pipeline does not apply a multiple-comparison correction. Compare models and biological contexts deliberately rather than ranking ratios alone.
The HTML opens with six sequencing panels: assigned read pairs, barcode balance within each physical library, mapped reads, usable fragments after filtering, PCR duplication, and end-to-end usable yield. The five peak panels show peak count, FRiP, total bases covered, the peak-width median and range, and peak count versus usable fragments. Values are printed above bars, and axes include their units.
Insert-size and peak-width distributions use normalized, shared 250-bp bins
whose percentages are comparable across samples. The underlying
insert_size_distribution and width_distribution arrays in
qc_summary.json remain unbinned. The fragments-per-peak coverage ECDF is
built from
qc/peaks/<sample_id>/<sample_id>.peak_qc.fragments_per_peak.tsv; its public
JSON representation is a compact fragment_count/peak_count histogram.
Peaks with no overlapping fragments are retained as zero-fragment peaks and
reported explicitly rather than dropped or converted to missing values.
Across bar, distribution, fragments-per-peak, TSS, and scatter panels, the target palette uses stable colors for IgG, CTCF, GATA1, and RUNX1, with deterministic fallback colors for other targets. Line plots use direct endpoint labels and distinct markers/dashes so the samples remain identifiable without color alone.
The motif enrichment heatmap consumes complete AME output. It displays all
cognate motifs detected by whole-token matching (for example, GATA1::TAL1
and TAL1::GATA1) plus the 15 non-cognate motifs with the strongest adjusted
p-values across the cohort. GATA10 therefore does not match an expected
GATA1. The separate top_motifs.tsv export remains capped at ten motifs per
target.
HTML tables and chart labels use three significant digits for compact presentation. Machine-readable TSV and JSON outputs retain full precision, so downstream calculations should consume those files rather than values copied from the HTML. Long supporting tables are collapsed by default.
The dashboard machine-readable contract is schema version 1 and generator version 1.1.0.
This is the first released form of the contract. Consumers
should check both values before interpreting fields. The HTML is a descriptive
view of these products, not a machine interface.
qc_summary.tsv has one row per derived sample, ordered by sample_id, with
these stable columns in this exact order:
sample_id, library_id, input_group, assay_target, is_control, control_id,
expected_motif, total_read_pairs, assigned_read_pairs, ambiguous_read_pairs,
unassigned_read_pairs, assigned_fraction, ambiguous_fraction,
unassigned_fraction, sample_assigned_reads, sample_assignment_fraction,
raw_total_reads, mapped_percent, properly_paired_percent, mapq_filtered_reads,
mapq_filtered_fragments, mapq_filtered_fraction, markdup_examined_reads,
duplicate_total, duplicate_percent, mitochondrial_percent,
estimated_library_size, insert_size_total_pairs, insert_size_min,
insert_size_q25, insert_size_mean, insert_size_median, insert_size_q75,
insert_size_max, peak_count, total_covered_bases, total_fragments,
fragments_in_peaks, frip, peak_width_min, peak_width_q25, peak_width_mean,
peak_width_median, peak_width_q75, peak_width_max, tss_status,
tss_enrichment, expected_motif_status, best_motif_id,
best_adjusted_p_value, ame_status, warning_count
top_motifs.tsv contains at most ten rank-ordered AME rows per target and
never contains IgG-control rows. Its stable columns are:
sample_id, assay_target, expected_motif, rank, motif_id, motif_alt_id,
adjusted_p_value, p_value, effect, positive_sequences
tss_profiles.tsv is the stable tidy profile export with columns
sample_id, position_bp, signal. The native
qc/tss/<sample_id>/<sample_id>.tss_profile.tsv is instead the pinned
deepTools 3.5.5 plotProfile --outFileNameData table and should not be treated
as the dashboard's stable consumer schema.
qc_summary.json has exactly these top-level fields:
schema_version, generator_version, annotation_status, counts,
availability, metric_definitions, samples, and warnings. counts
contains sample, target, control, and warning counts. availability summarizes
the demultiplex, library, insert_size, peak, peak_width, tss,
motif, and ame families. Each sample carries the same families with a
status and reason; statuses are computed, skipped, empty, missing,
failed, or not_applicable. Insert-size and peak-width distributions are
full-resolution, unbinned arrays under the corresponding sample's library
and peak objects. The target-only fragments_per_peak_distribution is a
compact histogram under peak. Producer-specific columns are not added
automatically to this versioned public object.
TSV missing numeric values are empty fields, booleans are lowercase true or
false, and finite numbers use locale-independent text. JSON missing values
are null; an observed but header-only distribution is an empty array with an
empty availability status. The HTML displays missing values as NA.
Absence, intentional skips, empty analyses, and failed/partial artifacts
therefore remain distinguishable from measured zero.
The QC subworkflow now requires eleven inputs, in order:
filtered_bams, final_broad_peaks, coverage, library_metrics,
demultiplex_metrics, fastqc_reports, motif_metrics, gtf, tss_bed,
motif_ame_results, and motif_ame_statuses. Direct importers of the previous
nine-input subworkflow must provide the two AME channels; use empty channels
when motif analysis is disabled.
The TSS_ENRICHMENT.out.profiles tuple now contains seven values:
meta, TSS BED, compressed matrix, matrix table, profile image, native
deepTools profile table, and status table. Direct module consumers that
destructure the earlier six-value tuple must insert the profile-table element
before the status table.
- Demultiplexing:
assigned_fraction,ambiguous_fraction, andunassigned_fractionuse all synchronized input read pairs as denominator. Ambiguous means two or more expected barcodes tied at the closest permitted Hamming distance. Unassigned means none passed the threshold. - Alignment and pairing:
mapped_percentandproperly_paired_percentare the percentages reported bysamtools flagstaton the primary analysis BAM. - MAPQ-filtered fraction: the filtered BAM retains properly paired primary
alignments (
samtools view -f 2 -F 2304) at or above--min_mapq.mapq_filtered_fractionis filtered alignment records divided by SAMtoolsraw total sequences; two validated mate records constitute one filtered fragment. - Duplicate fraction:
duplicate_percentis SAMtools markdupDUPLICATE TOTALdivided byEXAMINED(orREADfor compatible output), reported as a percentage.estimated_library_sizeis retained. Duplicates are measured, not removed from the primary analysis BAM. MACS2 intentionally uses--keep-dup 1; coverage independently uses--ignoreDuplicates. - Mitochondrial fraction: mapped records on
chrM,MT, orMdivided by mapped records on all named reference sequences, reported as a percentage. - FRiP: for each non-control sample, the numerator is the number of unique,
properly paired fragments from the filtered BAM overlapping at least one
final (blacklist-filtered when applicable) broad peak. The denominator is
all unique properly paired fragments in that filtered BAM. A fragment is
counted once even if both mates or multiple peaks overlap. The original FRiP
value remains at
qc/peaks/<sample_id>/<sample_id>.peak_qc.tsv, in rowfrip; the dashboard summary joins that source with the other QC products. - Peak QC: includes count, union-covered bases, width min/mean/median/max and quartiles, width histogram, MACS2 score/signal summaries, and fragment-per-peak counts. These describe the final broad peaks.
- TSS enrichment: when
--tss_bedor--gtfis supplied, deepTools builds a strand-aware matrix from 3 kb upstream to 3 kb downstream in 10-bp bins and emits the matrix and aggregate profile. The dashboard's scalar TSS enrichment is position 0 divided by the mean first/last 100 bp, and the full curve is retained intss_profiles.tsv. Without annotation it recordsskipped_no_annotation.
IgG controls are visually separated from targets in the dashboard. Controls have no expected-motif result, because expected motifs apply only to target libraries. When optional outputs are unavailable, the dashboard presents NA warnings rather than zero so absence is not mistaken for a measured value. The dashboard is descriptive and applies no biological thresholds; establish project-specific interpretation before making biological classifications.
These metrics are descriptive QC, not universal pass/fail thresholds. Compare targets with matched controls and comparable input groups, and establish project-specific cutoffs before excluding libraries.
The primary result reproduces the NanoScope-style, matched-IgG MACS2 BAMPE
call with --llocal 100000 --keep-dup 1 --broad-cutoff 0.1 --max-gap 1000 --broad. CUT&Tag signal can occupy extended regulatory domains, and preserving
this broad call is the compatibility objective. Final broad peaks drive
blacklist filtering, FRiP, peak QC, and the principal results.
Narrow calls use the same target/control BAMs, BAMPE mode, local lambda, and
duplicate setting, but exist only to provide localized summits for motif QC.
They never replace the broad results or FRiP. Set
--motif_use_narrow_peaks false to skip them and center motif windows on
broad-peak midpoints.
--motif_db must be a readable MEME-format DNA motif database. A vertebrate
JASPAR release converted to MEME format, such as
JASPAR2026_CORE_vertebrates.meme, is an appropriate starting point provided
its motif names/alternate names match the manifest regular expressions.
For each non-control sample, the pipeline:
- extracts reference-bounded windows around narrow summits by default, or final broad-peak midpoints;
- creates seeded, non-overlapping, approximately GC-matched genomic background that excludes final peaks and an optional blacklist;
- runs AME Fisher known-motif enrichment, STREME de novo discovery, and FIMO scanning;
- matches the manifest
expected_motifregular expression and reports motif rank, enrichment/effect statistic, significance value, hit count, fraction of peaks with a hit, hit-position distribution, and central-hit fraction.
Expected expressions for this dataset are CTCF, GATA1, and RUNX1.
Matching is a regular-expression search against motif ID and alternate ID, so
inspect the chosen database naming before a production run. motif_not_found
is an explicit QC failure when a complete AME database report contains no
match; not_significant and no_peaks are distinct statuses.
For modern AME E-value-only output, the compatibility field named
best_adjusted_p_value carries the AME E-value. Interpret it according to the
AME report rather than assuming the field always contains an adjusted
probability.
Nextflow stores task state in the work directory and .nextflow history. Keep
those paths and rerun the identical command with -resume:
nextflow run main.nf \
-profile conda \
-work-dir /shared/nanocut/work \
--input /shared/project/samples.csv \
--fasta /shared/references/GRCh38.fa \
--macs_genome_size hs \
--outdir /shared/project/results \
-resumeCompleted tasks are eligible for cache reuse only when their inputs, command,
configuration, and relevant parameters are unchanged. Moving/changing input
files, deleting the work directory, changing the work path, or changing
parameters can invalidate cache entries. NXF_CONDA_CACHEDIR reuses software
environments but is separate from the task cache. Do not delete work/ until
the run is accepted and no resume is needed.
- Manifest fails immediately: run
python3 bin/manifest.py validate --input samples.csv --output normalized.jsonfor the precise row-level error. Check the exact header, relative FASTQ paths, identifier characters, boolean spelling, repeated physical-library paths, and same-input IgG links. - FASTQs are unsynchronized: ensure R1, R2, and I2 contain the same records
in the same order.
/1and/2suffixes are normalized, but different core read IDs or truncated records fail the affected library. - Many reads are ambiguous/unassigned: confirm barcodes are in I2 FASTQ
orientation and compare observed-I2 counts in demultiplex metrics. Increase
--barcode_mismatchesonly deliberately; it must remain below barcode length, and tied closest barcodes remain ambiguous. - Zero-read sample: correct the manifest/barcode or use
--allow_empty trueonly for diagnostics. The default failure prevents silent downstream empty analyses. - Index error: pass the prefix without
.1.bt2; verify all six small or all six large index files exist and that both sets are not present under the same prefix. - Motifs do not run: supply both
--fastaand a valid--motif_db. Confirm the MEME header andMOTIFrecords, expected-name regex, usable peaks, and background-generation logs. - No TSS profile: supply a valid BED6
--tss_bedor a GTF with strand-bearingtranscriptfeatures. When both are supplied, BED wins. - Conda environment creation is repeated or fails on a cluster: ensure
NXF_CONDA_CACHEDIRis the same absolute shared path on all nodes, writable for creation, and not on node-local storage. Solve eachenvs/*.ymlindependently to diagnose package/channel issues. - Docker image failure: verify Docker access, Linux image architecture, registry/network access, and the process-specific image tag. Image pulls were not verified here.
- A resume reruns tasks: use the same
-work-dir, preserve.nextflowhistory, and compare inputs/parameters/configuration. Reviewpipeline_info/trace.txt.
Available checks are:
PYTHONDONTWRITEBYTECODE=1 python3 -m unittest discover \
-s tests/unit -p 'test_*.py' -q
for test_script in tests/integration/*.sh; do
bash -n "$test_script"
bash "$test_script"
doneIf pytest is installed, the unit suite is also conventionally discoverable
with pytest tests/unit -q. Runtime-aware integration scripts execute their
Nextflow harnesses when Nextflow is available and otherwise print an explicit
SKIP.
Current limitations:
- Nextflow runtime execution was not verified in this workspace, so no
successful DSL2, end-to-end, or
-resumeruntime claim is made. - Conda/Mamba/Micromamba were unavailable, so the pinned environments were structurally checked but not solved.
- Real Bowtie2, SAMtools, deepTools, MACS2, BEDTools, MEME Suite, and MultiQC execution was not performed here; direct fixtures and offline fakes do not replace a production-scale validation.
- Container tags are pinned but were not pulled or architecture-tested.
- The bundled
testprofile uses deterministic fake tools and tiny synthetic data; it verifies orchestration contracts, not biological correctness. - The workflow is specialized to the manifest assay targets
IgG,CTCF,GATA1, andRUNX1.