The development of PyPGx was heavily inspired by Stargazer __, another star-allele calling tool developed by Steven when he was in his PhD program at the University of Washington. Therefore, please also cite the following articles:

• Lee et al., 2019. Calling star alleles with Stargazer in 28 pharmacogenes with whole genome sequences __. Clinical Pharmacology & Therapeutics.
• Lee et al., 2018. Stargazer: a software tool for calling star alleles from next-generation sequencing data using CYP2D6 as a model __. Genetics in Medicine.

• Wroblewski et al., 2022. Pharmacogenetic variation in Neanderthals and Denisovans and implications for human health and response to medications __. bioRxiv.
• Botton et al., 2020. Phased Haplotype Resolution of the SLC6A4 Promoter Using Long-Read Single Molecule Real-Time (SMRT) Sequencing __. Genes.

If you find my work useful, please consider becoming a sponsor __.

* - Package - Anaconda - PyPI * - fuc - ✅ - ✅ * - scikit-learn - ✅ - ✅ * - openjdk - ✅ - ❌

There are various ways you can install PyPGx. The recommended way is via conda (Anaconda __):

Above will automatically download and install all the dependencies as well. Alternatively, you can use pip (PyPI __) to install PyPGx and all of its dependencies except openjdk (i.e. Java JDK must be installed separately):

The nice thing about this approach is that you will have access to development versions that are not available in Anaconda or PyPI. For example, you can access a development branch with the git checkout command. When you do this, please make sure your environment already has all the dependencies installed.

.. note:: Beagle __ is one of the default software tools used by PyPGx for haplotype phasing SNVs and indels. The program is freely available and published under the GNU General Public License license>_. Users do not need to download Beagle separately because a copy of the software (beagle.22Jul22.46e.jar) is already included in PyPGx.

Starting with the 0.12.0 version, reference haplotype panel files and structural variant classifier files in PyPGx are moved to the pypgx-bundle repository __ (only those files are moved; other files such as allele-table.csv and variant-table.csv are intact). Therefore, the user must clone the pypgx-bundle repository with matching PyPGx version to their home directory in order for PyPGx to correctly access the moved files (i.e. replace x.x.x with the version number of PyPGx you're using, such as 0.18.0):

Starting with version 0.22.0, you can now specify a custom location for the pypgx-bundle directory instead of using the home directory. This can be achieved by setting the bundle location using the PYPGX_BUNDLE environment variable:

Many pharmacogenes are known to have structural variation (SV) __ such as gene deletions, duplications, and hybrids. You can visit the Genes __ page to see the list of genes with SV.

Some of the SV events can be quite challenging to detect accurately with NGS data due to misalignment of sequence reads caused by sequence homology with other gene family members (e.g. CYP2D6 and CYP2D7). PyPGx attempts to address this issue by training a support vector machine (SVM) __-based multiclass classifier using the one-vs-rest strategy __ for each gene for each GRCh build. Each classifier is trained using copy number profiles of real NGS samples as well as simulated ones, including those from 1KGP __ and GeT-RM __.

PyPGx was recently applied to the entire high-coverage WGS dataset from 1KGP (N=2,504). Click here __ to see individual SV calls, and corresponding copy number profiles and allele fraction profiles.

When working with PGx data, it's not uncommon to encounter a situation where you are handling GRCh37 data in one project but GRCh38 in another. You may be tempted to use tools like LiftOver to convert GRCh37 to GRCh38, or vice versa, but deep down you know it's going to be a mess (and please don't do this). The good news is, PyPGx supports both of the builds!

In many PyPGx actions, you can simply indicate which genome build to use. For example, for GRCh38 data you can use --assembly GRCh38 in CLI and assembly='GRCh38' in API. **Note that GRCh37 will always be the default.** Below is an example of using the API:

However, there is one important caveat to consider if your sequencing data is GRCh38. That is, sequence reads must be aligned only to the main contigs (i.e. chr1, chr2, ..., chrX, chrY), and not to the alternative (ALT) contigs such as chr1KI270762v1alt. This is because the presence of ALT contigs reduces the sensitivity of variant calling and many other analyses including SV detection. Therefore, if you have sequencing data in GRCh38, make sure it's aligned to the main contigs only.

The only exception to above rule is the GSTT1 gene, which is located on chr22 for GRCh37 but on chr22KI270879v1alt for GRCh38. This gene is known to have an extremely high rate of gene deletion polymorphism in the population and thus requires SV analysis. Therefore, if you are interested in genotyping this gene with GRCh38 data, then you must include that contig when performing read alignment. To this end, you can easily filter your reference FASTA file before read alignment so that it only contains the main contigs plus the ALT contig. If you don't know how to do this, here's one way using the fuc program (which should have already been installed along with PyPGx):

In order to efficiently store and transfer data, PyPGx uses the ZIP archive file format (.zip) which supports lossless data compression. Each archive file created by PyPGx has a metadata file (metadata.txt) and a data file (e.g. data.tsv, data.vcf). A metadata file contains important information about the data file within the same archive, which is expressed as pairs of =-separated keys and values (e.g. Assembly=GRCh37):

* - Metadata - Description - Examples * - Assembly - Reference genome assembly. - GRCh37, GRCh38 * - Control - Control gene. - VDR, chr1:10000-20000 * - Gene - Target gene. - CYP2D6, GSTT1 * - Platform - Genotyping platform. - WGS, Targeted, Chip, LongRead * - Program - Name of the phasing program. - Beagle, SHAPEIT * - Samples - Samples used for inter-sample normalization. - NA07000,NA10854,NA11993 * - SemanticType - Semantic type of the archive. - CovFrame[CopyNumber], Model[CNV]

• CovFrame[CopyNumber]

• CovFrame[DepthOfCoverage]

• CovFrame[ReadDepth]

• Model[CNV]

• SampleTable[Alleles]

• SampleTable[CNVCalls]

• SampleTable[Genotypes]

• SampleTable[Phenotypes]

• SampleTable[Results]

• SampleTable[Statistics]

• VcfFrame[Consolidated]

• VcfFrame[Imported]

• VcfFrame[Phased]

To demonstrate how easy it is to work with PyPGx archive files, below we will show some examples. First, download an archive to play with, which has SampleTable[Results] as semantic type:

We can unzip it to extract files inside (note that tmpcty4c_cr is the original folder name):

Please visit the Genes __ page to see the list of genes with a genotype-phenotype table and each of their prediction method.

To perform phenotype prediction with the API, you can use the pypgx.predict_phenotype method:

To perform phenotype prediction with the CLI, you can use the call-phenotypes command. It takes a SampleTable[Genotypes] file as input and outputs a SampleTable[Phenotypes] file:

PyPGx currently provides three pipelines for performing PGx genotype analysis of single gene for one or multiple samples: NGS pipeline, chip pipeline, and long-read pipeline. In additional to genotyping, each pipeline will perform phenotype prediction based on genotype results. All pipelines are compatible with both GRCh37 and GRCh38 (e.g. for GRCh38 use --assembly GRCh38 in CLI and assembly='GRCh38' in API).

Implemented as pypgx run-ngs-pipeline in CLI and pypgx.pipeline.runngspipeline in API, this pipeline is designed for processing short-read data (e.g. Illumina). Users must specify whether the input data is from whole genome sequencing (WGS) or targeted sequencing (custom targeted panel sequencing or whole exome sequencing).

This pipeline supports SV detection based on copy number analysis for genes that are known to have SV. Therefore, if the target gene is associated with SV (e.g. CYP2D6) it's strongly recommended to provide a CovFrame[DepthOfCoverage] file and a SampleTable[Statistics] file in addtion to a VCF file containing SNVs/indels. If the target gene is not associated with SV (e.g. CYP3A5) providing a VCF file alone is enough. You can visit the Genes __ page to see the full list of genes with SV. For details on SV detection algorithm, please see the Structural variation detection __ section.

When creating a VCF file (containing SNVs/indels) from BAM files, users have a choice to either use the pypgx create-input-vcf command (strongly recommended) or a variant caller of their choice (e.g. GATK4 HaplotypeCaller). See the Variant caller choice __ section for detailed discussion on when to use either option.

Check out the GeT-RM WGS tutorial __ to see this pipeline in action.

Implemented as pypgx run-chip-pipeline in CLI and pypgx.pipeline.runchippipeline in API, this pipeline is designed for DNA chip data (e.g. Global Screening Array from Illumina). It's recommended to perform variant imputation on the input VCF prior to feeding it to the pipeline using a large reference haplotype panel (e.g. TOPMed Imputation Server __). Alternatively, it's possible to perform variant imputation with the 1000 Genomes Project (1KGP) data as reference within PyPGx using --impute in CLI and impute=True in API.

Check out the Coriell Affy tutorial __ to see this pipeline in action.

Implemented as pypgx run-long-read-pipeline in CLI and pypgx.pipeline.runlongread_pipeline in API, this pipeline is designed for long-read data (e.g. Pacific Biosciences and Oxford Nanopore Technologies). The input VCF must be phased using a read-backed haplotype phasing tool such as WhatsHap __.

PyPGx outputs per-sample genotype results in a table, which is stored in an archive file with the semantic type SampleTable[Results]. Below, we will use the CYP2D6 gene with GRCh37 as an example to illustrate how to interpret genotype results from PyPGx.

• Genotype: Diplotype call. When there is no SV this simply combines the two top-ranked star alleles from Haplotype1 and Haplotype2 with the delimiter '/'. In the presence of SV the final diplotype is determined using one of the genotypers in the pypgx.api.genotype module (e.g. CYP2D6Genotyper __).
• Phenotype: Phenotype call.
• Haplotype1, Haplotype2: List of candidate star alleles for each haplotype. For example, if a given haplotype contains three variants 22-42523943-A-G, 22-42524947-C-T, and 22-42526694-G-A, then it will get assigned 4;10; because the haplotype pattern can fit both \4 (22-42524947-C-T) and \10 (22-42523943-A-G and 22-42526694-G-A). Note that \4 comes first before \10 because it has higher priority for reporting purposes (see the pypgx.sortalleles method alleles>__ for detailed implementation).
• AlternativePhase: List of star alleles that could be missed due to potentially incorrect statistical phasing. For example, let's assume that statistical phasing has put 22-42526694-G-A for Haplotype1 and 22-42523805-C-T for Haplotype2. Even though the two variants are in trans orientation, PyPGx will also consider alternative phase in case the two variants are actually in cis orientation, resulting in 69; as AlternativePhase because \69 is defined by 22-42526694-G-A and 22-42523805-C-T.
• VariantData: Information for SNVs/indels used to define observed star alleles, including allele fraction which is important for allelic decomposition after identifying CNV (e.g. the sample NA19207). In some situations, there will not be any variants for a given star allele because the allele itself is "default" allele for the selected reference assembly (e.g. GRCh37 has \2 as default while GRCh38 has \1).
• CNV: Structural variation call. See the Structural variation detection __ section for more details.