audit examples/ dir -- what do we need?

[andrewkern]

Review the examples/ directory and determine what's needed for a complete, useful set of examples.

• Inventory current scripts and notebooks
• Identify gaps in coverage (e.g. missing use cases, modules without examples)
• Decide which existing benchmark scripts should stay vs. be moved/removed
• Ensure examples are minimal, well-documented, and serve as clear templates for users

[andrewkern]

Did an audit of what's in examples/. Here's the proposal:

Keep (5 files):

• moments_integration_demo.py — flagship demo, well-documented
• pg_gpu_tour.ipynb — interactive walkthrough, referenced in README
• genome_scan.py — practical real-world workflow
• multipop_ld_genotype_benchmark.py — covers 2/3/4-pop LD benchmarks
• vcf_to_zarr.py — useful utility

Remove (4 files):

• two_pop_ld_benchmark.py — superseded by multipop_ld_benchmark.py
• four_pop_ld_demo.py — superseded by multipop_ld_benchmark.py
• multipop_ld_benchmark.py — superseded by multipop_ld_genotype_benchmark.py (genotype version is more general)
• windowed_benchmark.py — benchmarks against scikit-allel which isn't the relevant comparison target anymore

Also considering moving the utility scripts (vcf_to_zarr.py, genome_scan.py) into a utils/ directory to separate them from demos and benchmarks.

Started work in #61.

[andrewkern]

also -- we could collapse multipop_ld_genotype_benchmark.py and multipop_ld_benchmark.py into one script and add a command line option that toggles between haploid and diploid modes

[andrewkern]

@nspope what do you think?

[nspope]

I largely agree with your suggestions. We don't need so many moments LD scripts--- really just a single benchmark and the integration demo. (You could keep the benchmark for various number of populations, but honestly the 4pop statistics alone should be sufficient, as they include the lower-order stats).

As far as what's missing---the docs and the tour notebook do a good job of showing general usage. I think other "examples" should consist of vignettes of various useful analyses than can be accelerated by pg_gpu. Some thoughts:

• admixture analysis? e.g. scanning for introgression, or detecting populations that are introgressed with f3/D.
• block bootstrapping? we could add infrastructure to do this easily, or just spell out how to do it in an example
• local pairwise relatedness calculation to local PCA?
• calculating isolation-by-distance curves? this should be very fast given the pairwise stats
• LD partitioning into blocks? Again the idea would be to use the fast pairwise stats, and hook it up to some existing graph based approach for delimiting blocks

[nspope]

... one future idea (not for now, not really relevant to this thread) is to do ANGSD-style site-frequency likelihood computation from genotype likelihoods. This is parallel across sites, and can also be blocked across groups of individuals (combining blocks is a convolution that is linear in the size of the block).

[andrewkern]

okay these are good suggestions. I'll work on trying to pilot some of them and see what looks clean.

[nspope]

another thing we should clearly show in an example: how to use the accessible mask to filter sequence prior to analysis. For example, removing an exon annotation prior to moments or admixture analysis. This could be incorporated into one of the existing examples.

[andrewkern]

that's an easy one! i can knock that out

[andrewkern]

• LD partitioning into blocks? Again the idea would be to use the fast pairwise stats, and hook it up to some existing graph based approach for delimiting blocks

@nspope, this seems like a great fit, but do you have a graph approach in mind? scipy.cluster.hierarchy ?