This skill drives postdoc-grade phylodynamic inference -- time-scaled phylogenies, R_e estimation, structured-population analysis -- from dated pathogen genomes. It covers ML phylodynamics (TreeTime), full Bayesian (BEAST2 with BICEPS, BDSKY, MASCOT, ORC clock), Bayesian bacterial dating (BactDating after Gubbins), pandemic-scale phylogenetics (UShER + matUtils + RIPPLES), and recombination masking for bacteria (Gubbins, ClonalFrameML). The decision-grade content is everything that surrounds the tool calls: TempEst temporal-signal QC first, recombination masking before any bacterial clock inference, BDSKY origin-vs-rootHeight discipline, multi-chain convergence for any BEAST report, MASCOT over BEAST DTA for source attribution, BICEPS over BSP as the modern skyline default, and reconciliation between phylodynamic R_e and case-based R_t with explicit sampling-bias accounting.
Distinct from phylogenetics/divergence-dating (general-purpose dating beyond outbreaks) and from phylogenetics/bayesian-inference (BEAST mechanics beyond BDSKY). This skill owns the surveillance application -- R_e estimation, Ne(t), source attribution, recombination-aware bacterial phylodynamics.
conda install -c bioconda beast2 treetime iqtree gubbins clonalframeml snippy snp-dists usher
packagemanager -add BDSKY BEASTLabs feast ORC MASCOT BICEPS
pip install dendropy biopython baltic
R -e "install.packages(c('BactDating', 'bdskytools', 'coda', 'ape', 'ggplot2'), repos='https://cran.r-project.org')"
Tell the AI agent what phylodynamic question to answer:
"Sixty Mycobacterium tuberculosis clinical isolates from a hospital cluster investigation spanning 18 months. Snippy against H37Rv, snippy-core for the core alignment, Gubbins on
core.full.alnfor recombination masking (even though Mtb is largely clonal, document recombination check for cross-lineage analysis), IQ-TREE on the masked alignment with+ASC, TreeTime for time-scaling with--coalescent skyline --clock-filter 4. Inspectroot_to_tip_regression.pdfand report R^2; if R^2 < 0.3, flag insufficient temporal signal and do not proceed to BDSKY. Run date-randomisation as secondary diagnostic. Document the bundled Coll/Napier barcode version used by TB-Profiler for lineage assignment."
"Two hundred Klebsiella pneumoniae isolates from a 3-year regional carbapenemase surveillance project. Snippy + snippy-core, Gubbins to mask recombination on
core.full.aln(Klebsiella is recombining; mandatory), IQ-TREE on the masked alignment with+ASC, then BEAST2 with BDSKY in BEAUti template. Setorigin > rootHeightper Stadler 2013 (init to tMRCA + 0.1*tMRCA). Set sampling proportion per quarter (the actual fraction of clinical isolates sequenced). Become-uninfectious rate from epi knowledge (1/14 days). Run 4 independent BEAST chains from seeds 42, 17, 99, 7; combine via logcombiner only after marginal posteriors overlap; report Gelman-Rubin via coda. Compare phylodynamic R_e to case-based R_t from local surveillance data."
"One thousand Streptococcus pneumoniae genomes from a national post-PCV13 carriage surveillance. Build core-genome alignment via snippy-core against ATCC 700669 reference; run Gubbins on
core.full.aln(S. pneumoniae recombines intensely; unmasked clock rate would inflate 2-5x); rebuild tree with IQ-TREE on the recombination-masked alignment; compare clock rate to Croucher 2013 literature (~1.5e-6 subs/site/year) as a sanity check. If rate is >5e-6, recombination masking is insufficient -- consider fastGEAR (Mostowy 2017) for the highly-recombinogenic lineages."
"Four hundred SARS-CoV-2 sequences from 5 European countries during Alpha emergence. >=20 sequences per country. BEAUti -> MASCOT template; one trait per country; if sampling intensity covaries with travel volume, switch to MASCOT-GLM. Document why MASCOT not BEAST DTA (Lemey 2009 DTA is sampling-biased per De Maio 2015). Report migration rate matrix with HPD intervals; flag any HPD spanning >2 orders of magnitude as identifiability failure."
"Place 5,000 new SARS-CoV-2 sequences on the UCSC daily-updated UShER MAT; use matUtils to extract the JN.1-descendant subtree; re-estimate branch lengths with TreeTime (UShER parsimony branch lengths are biased low); fit BDSKY on the subtree for lineage-specific R_e through time. Document that R_e is lineage-specific phylodynamic, NOT population-wide; reconcile with case-based R_t."
"These 80 dated influenza H3N2 HA sequences span 6 months of a single season. Run TempEst root-to-tip regression; if R^2 < 0.3, the data do not support time-scaled inference. Run date-randomisation (TreeTime
clock --reassign-dates) as secondary check; note that date-randomisation can pass with narrow sampling windows (false negative). If temporal signal insufficient, report a tree without time-scaling and explain why R_e estimation is not supported."
core.full.aln (NOT core.aln); rebuild tree on the recombination-masked alignment.--coalescent skyline --clock-filter 4 --confidence.origin > rootHeight for BDSKY; set sampling proportion per epoch.packagemanager; record random seed.run_gubbins.py core.aln is wrong. The input must be core.full.aln (full positions including invariant), or Gubbins cannot estimate background SNP density.origin is the time from the start of the epidemic to the most recent sample. rootHeight is the time back to the MRCA of the sample. origin > rootHeight always; setting them equal biases R_e upward systematically.