Skip to contents

Try popgenVCF in 5 minutes

Source: vignettes/quickstart.Rmd
quickstart.Rmd

This is the fastest way to see real popgenVCF output: no Docker, no download, no configuration file. The package ships a small but real dataset – 160 individuals from 8 populations in the 1000 Genomes Project’s chromosome 22, X, and Y callsets – and every analysis runs on it exactly as it would on your own VCF.

Prefer to look before running anything? Results and Interpretation shows real figures from this exact dataset, and the complete example PDF report is committed in the repository.

1. Locate the bundled dataset

See system.file("extdata", "quickstart", "README.md", package = "popgenVCF") for exactly which real samples and populations are included and why (it deliberately contains a known real duplicate/MZ-twin pair, so the kinship analysis below has genuine, verifiable signal – and, as a bonus, a cautionary real example of why a single flagged pair needs corroborating evidence, not blind acceptance; see “What to look for” below).

2. Run the pipeline 

cfg <- default_config()
cfg$input$vcf <- paths$vcf
cfg$input$metadata <- paths$metadata
cfg$output$directory <- tempfile("popgenvcf-quickstart-")
# cfg$compute$threads/memory_mb are already auto-detected with headroom
# reserved (one thread spared, memory capped at 85% of what's available) --
# no manual override needed.

analysis <- run_pipeline(cfg)

On a modern multi-core machine this takes a few minutes – most of the time is DAPC’s cross-validation across the default K = 2:10 range, run in full because it is exactly what a real analysis does. Every module runs with its documented defaults: QC and LD pruning, PCA, IBS, kinship, sex check, runs of homozygosity (with a length-class breakdown), FST (with its Nm gene-flow estimate, Jost’s D, and population-specific FST/beta), genome scans, linkage disequilibrium decay, LD-based effective population size, population genetic distance and tree, a population assignment test, diversity (including allelic richness and effective number of alleles), a site frequency spectrum with a mode-shift bottleneck screen, DAPC, AMOVA, and a neighbour-joining tree of individuals. Mantel/isolation-by-distance and the spatial autocorrelation correlogram run too: the bundled metadata includes real, documented population-level collection-site coordinates (see the dataset’s README.md for the exact source and caveats), so both analyses have real geographic distances to work with rather than being skipped.

3. Look at the output 

list.files(cfg$output$directory)
browseURL(file.path(cfg$output$directory, "report", "population_genomics_report.html"))

report/population_genomics_report.pdf is the compact, portable counterpart. Both are self-contained – copy either one without the rest of the result directory. tables/ and figures/ hold every individual result; analysis_results.rds is the complete serialized object for further work in R (get_analysis_result(analysis, "kinship"), and so on).

4. What to look for

  • Kinship (tables/36_close_relatives.tsv, figures/21_kinship_heatmap): both seeded high-kinship pairs show up classified duplicate/MZ twin with kinship near 0.5, but only one of them actually is a duplicate. NA19331/NA19334 (both LWK) are confirmed consistent: real, matching male chromosome X and Y genotypes, no contradicting evidence. HG03873/HG03998 (labelled as two different populations, ITU and STU) is the more interesting case: sex check below shows their real chromosome X and Y genotypes are genuinely different sexes, which rules out MZ twins/duplicate despite the high autosomal kinship – a real example of why a single flagged pair needs corroborating evidence, not blind acceptance.
  • Sex check (tables/42_sex_check.tsv, figures/27_sex_check_F_by_sample): genetic sex inferred from chromosome X heterozygosity and chromosome Y call rate together, compared against the bundled sex metadata column. When the two signals confidently disagree the combined call is discordant, not silently resolved – this is what actually reveals the HG03873/HG03998 kinship caveat above.
  • PCA (figures/07_PCA_PC1_PC2): the 8 populations should separate broadly by continental origin. figures/07b_PCA_PC1_PC2_by_sex is the same PCA fit recoloured by the bundled sex metadata column instead – one static panel per metadata column that qualifies (enough samples per value; see the interpretation guide for the exact rule). Sex should show no separation here, which is itself the expected, useful result: population structure, not sex, drives this PCA.
  • Runs of homozygosity (tables/38_ROH_sample_summary.tsv): FROH is scaled to this run’s ~1Mb analyzed footprint, not the whole genome – see the interpretation guide before comparing it to anything else. The same table also splits FROH by run-length class (short/intermediate/long, Ceballos et al. 2018); at this footprint size no run can reach the “long” class, a mechanical consequence of the ~1Mb region rather than a biological finding – see the interpretation guide.
  • Mantel/isolation by distance (tables/25_Mantel_IBD_summary.tsv, figures/12_isolation_by_distance): a real, statistically significant positive signal (Mantel r = 0.313, p = 0.001) – genetic distance increases with geographic distance across these 8 continentally diverse populations. A partial Mantel test controlling for population identity (r = 0.163, p = 0.001) shows the signal weakens but survives once population membership is accounted for – see the interpretation guide for what this means.
  • Spatial autocorrelation (tables/50_spatial_autocorrelation.tsv, figures/50_spatial_autocorrelation): a distance-class correlogram complementing the Mantel test above. Read with real caution here: the bundled coordinates are population-level representative points, not per-individual GPS, so the shortest distance classes are empty – see the interpretation guide for why.
  • LD-based effective population size (tables/45_Ne_LD.tsv): needs markers on at least two chromosomes to work at all, and this dataset’s retained autosomal marker set is chr22-only by design, so every population correctly reports an honest “insufficient data” skip here – see the interpretation guide for what a real multi-chromosome VCF produces instead.
  • Population genetic distance and tree (tables/46_population_genetic_distance.tsv, trees/population_Nei_neighbor_joining.nwk): Nei’s standard genetic distance between populations, unlike the sample-level kinship/tree results above. LWK/YRI (both African) are the closest pair; the largest distances are all African-vs-non-African pairs – real continental population structure.
  • Population assignment test (tables/47_population_assignment.tsv, figures/47_population_assignment): 66% of samples correctly self-assign to their recorded population from genotype alone; most of the rest are ITU/STU or LWK/YRI cross-assignments – exactly the same two closely related population pairs the FST/Nm and Nei’s-distance results above already flagged, a real cross-method confirmation rather than a method weakness.
  • Site frequency spectrum and bottleneck screen (tables/49_bottleneck_mode_shift.tsv, figures/48_site_frequency_spectrum): all 8 populations show a mode-shifted spectrum with this run’s default qc.maf = 0.05 filter – but this is largely a filtering artifact, not a real signal; see the interpretation guide for the unfiltered comparison that shows only 3 of 8 are genuinely mode-shifted.

5. Next steps

Once the shape of the output is familiar, move to your own data with the Getting Started guide, which covers Docker/Apptainer execution, metadata preparation, and configuration in full.