What Pharmacogenomics Can Tell You From Your Own Genome Data

Pharmacogenomics (also written pharmacy genomics, or pharmacogenetics when it concerns single genes) is the mapping from your germline DNA to how your body absorbs, activates, metabolizes, and clears specific drugs. In practice it means converting a small number of pharmacogenes into diplotypes written in star allele notation (CYP2C19 *1/*17), translating those into a metabolizer phenotype, and looking up the gene-drug pairs that have published, evidence-graded guidance. If you have whole-genome sequencing you can produce most of this yourself in an afternoon with open-source tools. The result is informational until a clinician and a CLIA-confirmed test are in the loop, and there are three or four specific places where the pipeline will hand you a confident-looking answer that is wrong.

Diplotypes are the unit and variants are not

The mistake engineers make first is treating pharmacogenes like GWAS hits: grep the VCF for rs4244285, call yourself a poor metabolizer, move on. Pharmacogene function is defined at the haplotype level. A star allele is a named combination of variants on one chromosome copy, curated by PharmVar. The CYP2C19*2 allele is defined by rs4244285 (c.681G>A), but *35, *36, and *37 are separate haplotypes with their own functional assignments. A variant you call in isolation tells you nothing about which chromosome it sits on. Two loss-of-function variants in trans (one on each copy) and two in cis (both on the same copy, the other copy intact) are different phenotypes entirely.

So the pipeline is: variants → phased haplotypes → diplotype → CPIC function assignment → activity score → phenotype. Each arrow can fail independently.

These genes are worth calling, because they have CPIC guidelines at the level where the evidence is graded A.

• CYP2C19
• CYP2D6
• CYP2C9
• VKORC1
• DPYD
• TPMT
• NUDT15
• SLCO1B1
• CYP3A5
• UGT1A1
• CYP4F2
• G6PD
• IFNL3
• RYR1
• CACNA1S
• HLA-B/HLA-A for hypersensitivity alleles

That is the useful core. The broader literature on genetic variation in drug response covers far more genes, most of them at evidence levels that do not change anything you would do ¹.

The pipeline we would run

Start with a GRCh38-aligned BAM or CRAM and a joint-called or single-sample VCF. Two things matter before you run anything: build and normalization.

# Left-align and split multiallelics; PharmCAT is strict about this
bcftools norm -m- -c ws -f GRCh38.fa input.vcf.gz -Oz -o norm.vcf.gz
bcftools index norm.vcf.gz

Then PharmCAT, which is the reference implementation for VCF-to-phenotype:

pharmcat_pipeline norm.vcf.gz -o out/ -reporterJson

The preprocessor extracts the ~250 positions PharmCAT cares about, checks REF alleles against the reference genome, and complains loudly if your VCF is GRCh37. It will not silently lift over. The refusal is correct behavior: several pharmacogene positions sit in regions where the two builds disagree about the reference allele, and a silent liftover flips a call.

The single most dangerous flag is the one that treats missing positions as homozygous reference. If your VCF only contains variant sites (not a gVCF with reference blocks), PharmCAT cannot distinguish “reference at this position” from “no read coverage here.” Filling in reference by assumption will happily hand you *1/*1 for a gene where the region simply failed to sequence. Use a gVCF, or pass the BAM-derived callable regions and check depth yourself:

samtools depth -a -r chr10:94760000-94860000 sample.cram | awk '$3<10' | wc -l

If more than a handful of positions in the gene body are under 10x, do not trust the diplotype.

CYP2D6 is the hard one

CYP2D6 accounts for a large share of clinically relevant pharmacogene calls and is the gene short-read sequencing handles worst. It sits next to two pseudogenes, CYP2D7 and CYP2D8, with which it shares roughly 90 to 95 percent sequence identity. It has whole-gene deletions (*5), duplications and multiplications (*1x2, *2x3), and hybrid alleles where part of the gene is CYP2D7-derived (*13, *36, *68). Standard BWA-MEM plus GATK in this region produces mismapped reads, spurious heterozygous calls at paralogous positions, and copy number it cannot see at all.

PharmCAT does not call CYP2D6 from a VCF by default for this reason. Use a dedicated caller and feed the result in:

# Cyrius: CYP2D6 star alleles from WGS, handles CNV and CYP2D6-CYP2D7 hybrids
python star_caller.py --manifest manifest.txt --genome 38 \
    --prefix sample --outDir out/ --threads 8

# Aldy: multi-gene, works from BAM, reports copy number configurations
aldy genotype -p illumina -g CYP2D6 -o sample.aldy sample.bam

Run both. When Cyrius and Aldy disagree, the truth is usually a hybrid or a multiplication, and you should treat the call as unresolved rather than picking the answer you like. If the gene matters for something specific, long-read sequencing (PacBio HiFi or ONT) over the locus resolves it properly, because a single read spans the whole gene and phasing is free. That is the approach we would take when a CYP2D6 call is going to inform a conversation with a physician.

The failure modes that produce a wrong answer quietly

Array genotyping misses everything it was not designed to see. A SNP array with 20 pharmacogene probes calls *1 by default for any allele it has no probe for, and *1 means normal function. Rare and population-specific variants are the ones that get silently converted into reassurance. This is the strongest argument for sequencing over targeted genotyping in this domain ².

DPYD is the case where getting it wrong has the highest consequence. The four CPIC-actionable variants (including c.1905+1G>A / rs3918290 and the HapB3 haplotype marked by c.1129-5923C>G) are intronic or splice-region. Exome sequencing frequently does not cover them. If your data is exome, check coverage at those exact coordinates before believing a normal-function result.

G6PD is X-linked with hundreds of variants, and heterozygous females have an indeterminate phenotype that no genotype call resolves. HLA alleles need imputation or direct typing (HLA-B*57:01, HLA-B*15:02, HLA-A*31:01). PharmCAT will report them only if you supply an outside call from a tool like HLA-LA or Optitype. And structural variation in SLCO1B1 and UGT1A1 (the *28 TA repeat is a length polymorphism in a promoter) is poorly resolved by short reads.

Finally: a germline result is a static prior. The result does not measure current enzyme activity. Phenoconversion is real. Inflammation, liver function, and co-administered inhibitors and inducers can push a genotypic normal metabolizer into a poor metabolizer phenotype. Genotype tells you the ceiling. Serial biomarkers, transcript levels, and plasma protein measurements tell you where you are sitting under it today.

Reading the report without over-reading it

For each gene, PharmCAT’s HTML report gives you the diplotype, the phenotype, and the activity score where applicable. It also gives the CPIC recommendation text for each drug. Read the “uncalled” and “possibly missing” sections first, before the recommendations. Those sections are where the information about your data quality lives.

Two boundaries we hold to. First, none of this is a basis for starting, stopping, or changing a dose of anything. CPIC guidance is written for prescribers who already have a clinical indication in front of them, and applying it in the absence of one inverts the logic. Take the report to a physician or a pharmacist with pharmacogenomics training, and expect them to order a confirmatory CLIA test before acting. Second, running WGS on a healthy person surfaces findings well outside pharmacogenomics, some of which carry real psychological and insurance consequences. The evidence base for population-scale sequencing of asymptomatic adults is still thin ³. Decide in advance what you want returned.

Questions people also ask

What is pharmaceutical genomics? Same field as pharmacogenomics: the study of how inherited sequence variation changes drug absorption, metabolism, transport, target binding, and immune response, plus the tooling that turns a person’s sequence into actionable gene-drug pairs ¹.

What is the future of pharmacogenomics? Three shifts are underway. Long-read sequencing turning CYP2D6 and other structurally complex loci from a research problem into a routine call. Pre-emptive panels, where the genotype is in the record before the prescription is written rather than after a bad response. And combining germline calls with concurrent measurements (liver enzymes, inflammatory markers, plasma proteomics) to estimate present metabolic capacity rather than the genetic ceiling alone.

Can I just use my 23andMe raw data? For a rough look, yes, and PharmCAT accepts converted array data. Treat every normal-function call as unconfirmed, because an array reports reference at any position it does not probe, and that is indistinguishable from a true reference call in the output file ².

How many genes are worth calling? Around 15 to 20 have CPIC guidelines with graded evidence and defined phenotype mappings. Everything beyond that is currently research-grade.

Footnotes

1. Sushant M. Ahire, Prerana S. Pawar, Vijayraj N. Sonawane, et al. Personalized Medicine and Pharmacogenomics: Impact of Genetic variation on Drug Response and Individualized Therapy - A Review. Asian Journal of Pharmacy and Technology, 2026. https://doi.org/10.52711/2231-5713.2026.00013 ↩ ↩²

2. Michael J. Wagner. Pharmacogenetics and Personal Genomes. Personalized Medicine, 2009. https://doi.org/10.2217/pme.09.55 ↩ ↩²

3. Noralane M. Lindor, Stephen N. Thibodeau, Wylie Burke. Whole-Genome Sequencing in Healthy People. Mayo Clinic Proceedings, 2017. https://doi.org/10.1016/j.mayocp.2016.10.019 ↩