bedpull 0.3.0

bedpull - Pull the query sequence from bam or fasta references using a bed file
Documentation
#!/usr/bin/env bash
# bedpull benchmark configuration.
#
# Copy this file to config.env (gitignored) and fill in the paths below, then
# `source config.env` before running any script in this directory, or pass
# `--config config.env` to run_all.sh.
#
# Paths already known to be correct for this machine are filled in below;
# everything else is a placeholder you need to update.

# --- Reference genomes (already local per the user) ---
export HS1_FASTA="/path/to/chm13v2.0.fa"
export HG38_FASTA="/path/to/hg38.fa"

# --- HG002 T2T-Q100 v1.1 diploid assembly (query for PAF mode) ---
# Ships gzipped; bedpull/minimap2/samtools all want an uncompressed or
# bgzip-recompressed FASTA with a .fai index. Decompress once:
#   zcat hg002v1.1.fasta.gz > hg002v1.1.fasta && samtools faidx hg002v1.1.fasta
export HG002_ASSEMBLY_FASTA="/path/to/hg002/hg002v1.1.fasta"

# --- minimap2 PAF alignments: HG002 haplotype contigs -> HS1 ---
# Build with, e.g.:
#   samtools faidx "$HG002_ASSEMBLY_FASTA" $(grep '_PATERNAL' "$HG002_ASSEMBLY_FASTA.fai" | cut -f1) > hg002.paternal.fa
#   minimap2 -cx asm5 --cs=long -t 16 "$HS1_FASTA" hg002.paternal.fa > hg002pat_to_hs1.paf
# (same for _MATERNAL). Both are needed — see README.md's note on why we try
# both haplotypes per window instead of trusting VCF GT phase to assembly
# haplotype identity.
export HG002_PATERNAL_TO_HS1_PAF="/path/to/hg002pat_to_hs1.paf"
export HG002_MATERNAL_TO_HS1_PAF="/path/to/hg002mat_to_hs1.paf"

# --- Official UCSC chain file: hg38<->hs1 confirmation leg ---
# Download from https://hgdownload.soe.ucsc.edu/goldenPath/hg38/liftOver/
# UCSC only publishes this one direction (not hs1ToHg38), so the confirmation
# leg lifts hg38 -> hs1: windows are called directly from the hg38<->hs1 PAF
# alignment (call_alignment_svs.py), lifted onto hs1, sequence extracted from
# hs1. No building required for the chain itself.
export HG38_TO_HS1_CHAIN="/path/to/hg38ToHs1.over.chain.gz"

# PAF for the confirmation leg: target=hg38 (matches the windows' coordinate
# system and is what call_alignment_svs.py's truth is called from), query=hs1
# (what gets extracted) — see README.md's "hs1<->hg38 confirmation leg"
# section for the exact command.
export HS1_TO_HG38_PAF="/path/to/hs1_to_hg38.paf"

# --- Self-built chain: hs1<->HG002 leg (see build_chain_from_paf.sh) ---
export HG002_PATERNAL_TO_HS1_CHAIN="$PWD/out/hg002pat_to_hs1.chain"
export HG002_MATERNAL_TO_HS1_CHAIN="$PWD/out/hg002mat_to_hs1.chain"

# --- GIAB HG002 v5.0q structural variant truth set (main leg only — the
# confirmation leg's truth comes from call_alignment_svs.py instead, since
# neither hs1 nor hg38 is HG002's own assembly) ---
export GIAB_HS1_STVAR_VCF="/path/to/HG002_CHM13v2.0_v5.0q_stvar.vcf.gz"
export GIAB_HS1_STVAR_BED="/path/to/HG002_CHM13v2.0_v5.0q_stvar.benchmark.bed"

# --- Tools (override only if not on PATH) ---
export BCFTOOLS_BIN="${BCFTOOLS_BIN:-bcftools}"
export SAMTOOLS_BIN="${SAMTOOLS_BIN:-samtools}"
export MINIMAP2_BIN="${MINIMAP2_BIN:-minimap2}"
export LIFTOVER_BIN="${LIFTOVER_BIN:-liftOver}"
# Only needed for the optional paftools.js liftover extra-tool comparison —
# see README.md's "Optional third tool" section.
export K8_BIN="${K8_BIN:-k8}"
export PAFTOOLS_JS="${PAFTOOLS_JS:-paftools.js}"
# Only needed for the optional pslMap / liftOver-multiple extra-tool
# comparisons — see README.md's "Two more UCSC kent tools" section.
export BEDTOPSL_BIN="${BEDTOPSL_BIN:-bedToPsl}"
export PSLMAP_BIN="${PSLMAP_BIN:-pslMap}"
export LIFTOVERMERGE_BIN="${LIFTOVERMERGE_BIN:-liftOverMerge}"
export BEDPULL_BIN="${BEDPULL_BIN:-$PWD/../target/release/bedpull}"

# --- Window selection defaults (accuracy leg) ---
export SV_FLANK_BP=500
export SV_MIN_LEN=50
export SV_MAX_WINDOWS=0   # 0 = no cap, use every qualifying SV

# --- Output directory ---
export BENCH_OUTDIR="$PWD/out"