Open Source DNA Analysis: What You Can Run Yourself

Yes, the entire analysis stack is open source and free, and you can run all of it on a laptop or a single cloud instance. Alignment (bwa-mem2, minimap2), variant calling (GATK HaplotypeCaller, DeepVariant, bcftools), annotation (Ensembl VEP, snpEff, bcftools csq), and query layers (GEMINI, Hail, plink2) are all published, maintained, and free for personal use. What is not free is the data. If you upload a 23andMe or AncestryDNA text file to Promethease or Genetic Genie, the software is doing a competent job on an input that covers roughly 0.02% of your genome at a few hundred thousand pre-chosen positions. The limiting factor in personal genomics is almost never the analysis software. It is whether you have sequence data worth analyzing.

Array raw data versus sequence data

A consumer array file is a tab-delimited list of about 600,000 to 900,000 rsIDs with a genotype call each: rs4988235 2 136608646 AA. It is small, it loads in pandas in a second, and it supports exactly the analyses that array genotypes support: ancestry inference, relative matching, common-variant polygenic scores, and lookups of well-characterized common SNPs. It cannot tell you about a variant that is not on the chip, which is essentially every rare variant, and it cannot tell you about indels, structural variants, copy number, or repeat expansions.

The failure mode people hit hardest is the rare pathogenic variant lookup. Arrays are optimized for common variants where cluster separation is clean. For rare variants the clusters are sparse and a single mis-assigned point produces a heterozygous call in someone who carries nothing. This is why clinical genetics services do not accept array-derived pathogenic calls without orthogonal confirmation in a CLIA-certified laboratory, and why you should treat any “pathogenic” hit from a raw-data upload service as a hypothesis to bring to a genetic counselor rather than a result.

Sequence data changes the question you can ask. A 30x whole genome is roughly 100 GB of gzipped FASTQ, yields about 4 to 5 million variants against the reference, and covers the coding regions of all ~20,000 genes plus everything between them. The analytical framing shifts from looking up a list of positions to filtering a variant set, which is the framing the field has used since population-scale sequencing became routine ¹.

The pipeline we would run

Start from FASTQ if you have it, because it lets you choose your reference and your caller. Use the GRCh38 full analysis set with decoys and alt contigs (GRCh38_full_analysis_set_plus_decoy_hla.fa), not a bare primary assembly. The decoy sequences absorb reads from sequences absent in the reference and cut false positives in the regions that would otherwise collect them.

bwa-mem2 mem -t 32 -K 100000000 -Y \
  -R '@RG\tID:L1\tSM:sample\tPL:ILLUMINA\tLB:lib1' \
  GRCh38_full_analysis_set_plus_decoy_hla.fa r1.fq.gz r2.fq.gz \
| samtools sort -@ 8 -m 2G -o sample.bam -
samtools index -@ 8 sample.bam

The read group line matters: most downstream callers refuse to run without SM. Burrows-Wheeler indexing is what makes this tractable at all, reducing the memory footprint of a human-genome index to a few gigabytes and the alignment of a short read to milliseconds ².

For calling, we use DeepVariant rather than GATK HaplotypeCaller for a single sample. It has fewer knobs, needs no base quality score recalibration step, and its convolutional model handles indels in homopolymers better than a hand-tuned haplotype model.

run_deepvariant \
  --model_type=WGS \
  --ref=GRCh38_full_analysis_set_plus_decoy_hla.fa \
  --reads=sample.bam \
  --output_vcf=sample.vcf.gz \
  --output_gvcf=sample.g.vcf.gz \
  --num_shards=32

Keep the gVCF. The plain VCF cannot distinguish “reference at this position” from “no coverage here,” and that distinction is the difference between a negative result and no result. If you want to understand what the caller is doing underneath, the standard reference on the mechanics of read pileups, genotype likelihoods, and quality filtering is still worth reading before you trust your own output ³.

Then normalize, because annotation tools compare strings and unnormalized indels silently fail to match databases:

bcftools norm -m -any -f ref.fa -Oz -o sample.norm.vcf.gz sample.vcf.gz
bcftools index -t sample.norm.vcf.gz

Annotation is where the real work is

Ensembl VEP is the tool we recommend, run offline against a downloaded cache so that nothing leaves your machine:

vep --offline --cache --dir_cache ~/.vep --assembly GRCh38 \
  --fasta ref.fa --everything --pick_allele_gene \
  --plugin dbNSFP,/data/dbNSFP4.7a_grch38.gz,REVEL_score,CADD_phred,AlphaMissense_score \
  --custom /data/clinvar.vcf.gz,ClinVar,vcf,exact,0,CLNSIG,CLNREVSTAT \
  --custom /data/gnomad.genomes.v4.1.sites.vcf.bgz,gnomADg,vcf,exact,0,AF \
  -i sample.norm.vcf.gz -o sample.vep.vcf.gz --vcf --compress_output bgzip

--pick_allele_gene chooses one consequence per allele per gene and prevents the output from exploding across transcripts. The ClinVar CLNREVSTAT field is the one to filter on: a two-star or better review status means multiple submitters agree, and one-star or “no assertion criteria provided” entries include a great deal of material that has since been reclassified. Publicly curated annotation of this kind was proposed as the way genome interpretation would have to scale, and that prediction has held ⁴.

For querying, load the annotated VCF into a database rather than grepping it. GEMINI builds a SQLite store where variants, genotypes, and annotations are joined and you can write ordinary SQL against them ⁵:

gemini load -v sample.vep.vcf.gz -t VEP --cores 8 sample.db
gemini query -q "select chrom, start, ref, alt, gene, impact, clinvar_sig, aaf_gnomad_all \
  from variants where impact_severity='HIGH' and aaf_gnomad_all < 0.001" sample.db

A single 30x genome typically leaves you a few hundred rare high-impact variants after that filter, most of them in genes with no established phenotype, plus loss-of-function calls in olfactory receptors and other tolerant gene families. Filtering further by gene panel (ACMG secondary findings list, or a panel relevant to a documented family history) is the only way to get the list down to something interpretable, and interpretation past that point belongs with a clinical geneticist.

Specialist tools for regions generic callers miss

Short-read callers do badly in segmental duplications and highly polymorphic loci, and several of those loci are the ones people most want to see. Cyrius genotypes CYP2D6 by handling the CYP2D7 paralog explicitly. StellarPGx and PyPGx cover a wider set of pharmacogenes. ExpansionHunter calls repeat expansions at known loci from the BAM. Manta and Delly call structural variants. HLA-LA types HLA alleles from the alignment. Any pharmacogenomic call you generate is a research result until a clinical laboratory confirms it, and the decision about what, if anything, it means for medication is a conversation with a physician or pharmacist.

Failure modes that quietly produce wrong answers

Chromosome naming is the most common. GRCh38 files from Ensembl use 1, files from UCSC and GATK use chr1, and a mismatch makes bcftools annotate silently annotate nothing. Check with bcftools index -s before you blame the data. Build mismatch is next: a variant position from a 2012 paper is GRCh37, and lifting over with CrossMap or picard LiftoverVcf will drop or misplace a small fraction of sites, so record which build every coordinate came from. Array files add strand ambiguity, where A/T and C/G SNPs cannot be resolved by allele matching alone and require frequency-based checking. And imputation, whether through Beagle locally or a public imputation server, produces dosages with an r² per variant. Anything below about 0.8 should not be used for a single-variant claim about one person, though it is fine in aggregate for a polygenic score.

Questions people also ask

Is there a free DNA analyzer available? Yes. Ensembl VEP, snpEff, bcftools, DeepVariant, GEMINI, plink2, and IGV are all free and open source, and the annotation sources they depend on (ClinVar, gnomAD, dbSNP, the PGS Catalog) are free to download. Promethease charges a small fee and Genetic Genie is free, but both are built on array files and inherit that limitation.

How can I check my DNA for free? If you already have a consumer array file, you can load it in Python, convert it to VCF with plink2 --23file, and annotate it with VEP offline at no cost. If you do not have any data, the sequencing itself is the part you pay for, and no amount of open source tooling substitutes for it.

Can I do a DNA test on my phone? No. A phone can display results from a test run elsewhere, and there are apps that parse array files, but sample preparation and sequencing require laboratory instruments. Portable nanopore devices like the MinION are laptop-scale rather than phone-scale and still need a wet-lab library prep.

Should I use GRCh38 or the T2T-CHM13 assembly? Use GRCh38 as your primary, because nearly every annotation database is keyed to it. Run T2T-CHM13 as a secondary alignment if you care about centromeric, acrocentric, or segmental-duplication regions, which GRCh38 either lacks or represents poorly.

Oak builds longitudinal molecular profiles of individuals: whole-genome sequencing, RNA sequencing, proteomics, blood biomarkers, and continuous glucose data, integrated into one model of you. Build your profile.

Footnotes

1. Tuuli Lappalainen, Alexandra J. Scott, Margot Brandt, et al. Genomic Analysis in the Age of Human Genome Sequencing. Cell, 2019. https://doi.org/10.1016/j.cell.2019.02.032 ↩

2. Ben Langmead, Cole Trapnell, Mihai Pop, et al. Ultrafast and memory-efficient alignment of short DNA sequences to the human genome. Genome Biology, 2009. https://doi.org/10.1186/gb-2009-10-3-r25 ↩

3. André Altmann, Peter Weber, Daniel Bader, et al. A beginners guide to SNP calling from high-throughput DNA-sequencing data. Human Genetics, 2012. https://doi.org/10.1007/s00439-012-1213-z ↩

4. Tim Hubbard, Ewan Birney. Open annotation offers a democratic solution to genome sequencing. Nature, 2000. https://doi.org/10.1038/35002770 ↩

5. Umadevi Paila, Brad A. Chapman, Rory Kirchner, et al. GEMINI: Integrative Exploration of Genetic Variation and Genome Annotations. PLoS Computational Biology, 2013. https://doi.org/10.1371/journal.pcbi.1003153 ↩