Pharmacogenetic Associations with ADME Variants and Virologic Response to an Initial HAART Regimen in HIV-Infected Women

Background Clinical response to highly active antiretroviral therapy (HAART) varies among different populations. A portion of this variability may be due to variation in genes involved in the absorption, distribution, metabolism, and excretion (ADME) of HAART. Design To identify genetic factors involved in virologic responses to HAART, 13 genes in ADME pathways were analyzed in a cohort of HIV-infected women on HAART. A total of 569 HIV-positive participants from the Women’s Interagency HIV Study who initiated HAART from 1994–2012 and had genotype data were included in these analyses. Methods Admixture maximum likelihood burden testing was used to evaluate gene-level associations between common genetic variation and virologic response (achieving <80 viral copies/mL) to HAART overall and with specific drug classes. Results: Six statistically significant (P<0.05) gene-level burden tests were observed with response to specific regimen types. CYP2B6, CYP2C19 and CYP2C9 were significantly associated with response to protease inhibitor (PI)-based regimens. CYP2C9, ADH1A and UGT1A1 were significantly associated with response to triple nucleoside reverse transcriptase inhibitor (NRTI) treatment. Conclusions Although no genome-wide associations with virologic response to HAART overall were detected in this cohort of HIV-infected women, more statistically significant gene-level burden tests were observed than would be expected by chance (two and a half expected, six observed). It is likely that variation in one of the significant genes is associated with virologic response to certain HAART regimens. Further characterization of the genes associated with response to PI-based treatment is warranted.


Introduction
Morbidity and mortality has declined markedly with the advent and widespread use of highly active antiretroviral therapy (HAART) for individuals with HIV infection [1,2]. While the primary factor driving response to HAART is adherence, virologic response to HAART varies among adherent individuals and across populations [3]. Variation in genes that affect the absorption, distribution, metabolism, and excretion (ADME) of antiretrovirals (ARV) is thought to account for a portion of the variability in responses to treatment [4][5][6].
Several studies among HIV-infected individuals have evaluated virologic and immunologic responses, as well as adverse effects such as hypersensitivity reactions, neurotoxicity, hepatotoxicity, and hyperbilirubinemia, in relationship to genetic polymorphisms [4,5]. The genes which have been most frequently examined in pharmacogenetic association studies in HIV disease include those in the cytochrome p450 family (CYP), as well as the ABCB1 and UGT1A1 genes. The protein products of these genes are reported to be involved in the metabolism and transport of two ARV drug classes, non-nucleoside reverse transcriptase inhibitors (NNRTIs) and protease inhibitors (PIs) [6][7][8].
Numerous studies have evaluated CYP-based genetic determinants of ARV exposure, but whether these polymorphisms translate to drug efficacy is controversial [9][10][11][12][13][14][15][16][17][18][19]. There has also been wide interest in the potential genetic influence of ABCB1 on HAART responses owing to the crucial role of P-glycoprotein (encoded by ABCB1) in the distribution and excretion of PIs and NNRTIs [7,[20][21][22][23]. However, variation in virologic responses to PIcontaining regimens in HIV-infected patients who carry variant alleles of two independent ABCB1 polymorphisms has not been consistently observed [20][21][22][23]. Enhanced immunologic recovery in carriers of these ABCB1 variants has been observed in a recent analysis of HIV infected Chinese persons treated with PI-containing regimens (n=275) [24]. Loss of function alleles in the UGT1A1 gene are associated with drug-related hyperbilirubinemia; carriers of this variant treated with the PI atazanavir develop more severe hyperbilirubinemia, since atazanavir acts as an inhibitor to UGT1A1 [25].
Pregnane × receptor (NR1I2) is a nuclear receptor activated by several PIs, the nucleoside reverse transcriptase inhibitor (NRTI) abacavir, and the NNRTI efavirenz [26][27][28]. NR1I2 regulates the expression of CYP3A4, CYP2B6, and ABCB1, which mediate the metabolism and transport of several NNRTIs and PIs [26]. Pharmacokinetic associations between NR1I2 variation and the PI atazanavir have been reported [26][27][28]; however, pharmacogenetic associations with virologic response to HAART have not been explored. NRTIs are believed to be metabolized through alcohol dehydrogenase and the alcohol dehydrogenase 1A gene (ADH1A) has been examined for gene-environment interactions in the metabolism of abacavir and alcohol, although clinically significant interactions have not been identified [29].
Despite the potential for personalizing HAART regimens based on genetic predictors of response, pharmacogenetic studies of HAART responses to date have not provided conclusive evidence of associations. Previous pharmacogenetic studies in treated HIVinfected populations have focused on potential functional single nucelotide polymorphisms (SNPs) in genes that metabolize NNRTIs and PIs; the aim of this study was to comprehensively investigate genetic variation in ADME genes in the Women's Interagency HIV Study (WIHS), a prospective multi-site observational study of multi-ethnic women infected with HIV, in relationship to virologic response. Associations between polymorphisms, derived from genomewide genotyping data, in CYP1A2, CYP2A6, CYP2B6, CYP2C9, CYP2C19, CYP2D6, CYP2E1, CYP3A4, CYP3A5, ADH1A, ABCB1, UGT1A1, and NR1I2 and virologic response to an intial HAART regimen were explored.

Women's Interagency HIV Study (WIHS) Study Design
The WIHS is a prospective study of HIV-infected women and a comparison group of HIVuninfected women who were recruited in four phases from six sites across the U.S.: Bronx/ Manhattan, New York; Brooklyn, New York; Washington, D.C.; Los Angeles, CA; San Francisco/Bay Area, CA; and Chicago, IL [30,31]. Data in these analyses were restricted to the first two recruitment waves and include HIV-infected women only. The original recruitment phase was conducted in 1994-1995 and the second phase of recruitment was done in 2001-2002. The HIV-infected women in the WIHS are representative of the racial/ ethnic population of HIV-infected women in the U.S. [30,31]. A more detailed description of the WIHS cohort has been published [30,31].
Participants are seen for in-person visits every six months during which trained medical interviewers administer an extensive questionnaire, a clinical exam is performed, and biological samples are collected. Laboratory testing includes the measurement of HIV-1 viral loads (copies/ml) and CD4+ cell counts. At each biannual visit, detailed information on participants' current treatment regimen, the ARV medications taken in the previous six months and medications used for comorbidities is obtained. Commencing at visit 9 (1998), data were collected on self-reported adherence to ARV medications as determined by a visual analog scale (100%, 95-99%, 75-95%, <75% over the past six months).

Inclusion Criteria for Current Analysis
HIV-infected women who consented to genetic studies, initiated a regimen of three or more ARVs during study follow-up, maintained that regimen for a minimum of two consecutive study visits, had HIV RNA measurements at the visit immediately prior to, the visit of, and the visit immediately after initiating a three or more drug regimen were eligible for this pharmacogenetic study. Participants who were ARV therapy experienced prior to initiating a three drug regimen were included in the study only if at least two new ARV drugs were implemented as part of the initiated three drug HAART regimen. From October 1994 to September 2012, 1,307 WIHS participants initiated HAART. A total of 632 women were ineligible for this analysis: 138 had missing viral load data, 63 had undetectable viral loads at all relevant visits, 262 discontinued HAART during the study period, six had substantial missing covariate data, and for 163 participants only one new ARV was added to an existing regimen to achieve the definition of HAART. From a total of 675 eligible HAART initiators, 51 women were excluded because their treatment regimen did not fall within one of the three treatment categories (NRTI-, NNRTI-, and PI-based regimens) investigated here and 55 women were excluded due to missing genotype data, leaving a final analytic dataset of 569 WIHS participants.

Outcome
HIV viral load levels were quantified with a nucleic acid sequence based amplification assay with a lower limit of quantification of 80 copies/ml. A positive virologic response ("responders") was defined as achievement of an undetectable viral load at the visit during which the HAART regimen was first reported or at the visit subsequent to it, which corresponds to a maximum of 54 weeks of treatment. Participants not achieving an undetectable viral load at the visit of first reported HAART regimen use (corresponding to the previous six months) or the subsequent visit were classified as "non-responders".

Genotyping
A subset of 13 genes was selected from whole-genome genotyping data generated from a cohort-wide GWAS utilizing the Omni Human 2.5M SNP Array (Illumina, San Diego CA). Genotypes for SNPs in genes involved in ADME pathways (e.g. CYP1A2, CYP2A6, CYP2B6, CYP2C9, CYP2C19, CYP2D6, CYP2E1, CYP3A4, CYP3A5, ADH1A, ABCB1, UGT1A1, and NR1I2) were extracted for participants eligible for this study. Genotype imputation was performed using IMPUTE2, with imputed genotypes generated by comparing our GWAS SNP data for each gene region with reference haplotype data from the 1000 Genomes Phase 1 variant set release (version 3). 90 Markov chain Monte Carlo iterations were used in the imputation process. Only genotypes with an imputation r 2 of 0.50 or higher are included in the analysis. Coverage for each of the gene regions was calculated in each racial group by calculating pairwise r 2 between our genotyped and imputed (imputation r 2 ≥ 0.5 or greater, MAF ≥ 0.05) SNPs and every 1000 Genomes SNP in the region (MAF ≥ 0.05) using the 1000 Genomes samples of each respective racial group. SNPs that had an r 2 of 0.95 or greater with our genotyped and imputed SNPs were considered "covered," and the proportion of SNPs covered in each region was calculated.

Covariates
Genetic ancestry for the entire cohort was quantified by performing principal components analysis on ancestry informative markers included on the genome-wide panel. Genetic ancestry principal components (PC) covariates were used to control for population stratification in association models. Adherence data for this analysis were assessed at the visit at which the participant achieved an undetectable viral load (visit of HAART initiation or post-HAART initiation) since the adherence variable at this visit reflects HAART adherence in the six months preceding outcome ascertainment. Adherence data for nonresponders was taken at the study endpoint, which was defined as the visit immediately after the HAART initiation visit. For modeling purposes, adherence was dichotomized as ≥ 95% or < 95% adherent. Type of HAART regimen used by the participants was categorized as dual NRTI/PI, dual NRTI/NNRTI, or ≥ 3 NRTIs.

Statistical Analysis
The analysis includes women who self-identified as African-American, Hispanic, or Non-Hispanic White. Logistic regression was used to test associations between each SNP and achievement of an undetectable viral load. Genotypes were modeled as an ordinal variable, where common allele homozygotes, heterozygotes and minor allele homozygotes were coded as 0, 1, and 2, respectively.
Adherence, pre-HAART viral load and nadir CD4+ count prior to initiation were evaluated as potential confounders. Ultimately, since SNPs showed no association with self-reported adherence, pre-HAART viral load or nadir CD4+ cell counts, these variables were not included in the final model. Analyses were performed with adjustment for the top nine genetic ancestry PCs. Stratification on self-reported race/ethnicity did not affect the results.
There is evidence that ARV drug classes vary in their affinity for particular enzymes, and as such, the analyses were carried out with all regimens combined as well as separately by triple NRTI-alone regimens, NNRTI-containing or PI-containing regimens.
Given the lack of independence among our SNPs, an admixture maximum likelihood (AML) burden test was employed to evaluate whether a greater set of nominal associations per gene were observed than would be expected under the null distribution. Essentially, the AML method postulates that a given proportion of SNPs (α) within a set of SNPs is associated with the outcome and that the magnitude of each associated SNP's effect will fall on a noncentral χ2 distribution with non-centrality parameter η, a measure closely related to that SNP's contribution to the genetic variance of the outcome variable. The α and η parameters are estimated using a pseudo-maximum likelihood method [32]. This burden test was run using the AMLcalc program (http://ccge.medschl.cam.ac.uk/software/aml) with 1000 simulations; and the maximum proportion of associated SNPs was set to 0.2 on the genotyped and imputed data. We also adjusted for the first nine ancestry principal components. A total of 52 burden tests were run; comprising four HAART analytic groups (all regimens and each of the three regimen types) with each of the 13 genes. Table 1 summarizes the characteristics of HAART initiators from October, 1994 to September 2012 in the WIHS by virologic responder status. Of the 569 women included in the analysis, 59% were responders (n=335). Among the 313 (55%) women for whom selfreported adherence data were available, 12.1% of the responders reported <95% adherence compared to 28.3% of the non-responders.

Results
The 13 ADME genes evaluated in this analysis are shown in Table 2. The percentage of the genome covered in these regions by the genotyped and imputed SNPs by racial/ethnic group ranged from 17% to 100%. CYP2A6 and CYP3A5 were poorly covered in all of the racial/ ethnic groups. Although the genotyping array used in this project was designed to provide good coverage for non-Whites, only six of the 13 genes were covered at 90% or higher for African-Americans compared to 10 and nine genes, respectively for Whites and Hispanics.
When considering the three types of HAART regimens (triple NRTIs, dual NRTI/PI, dual NRTI/NNRTI), individual associations at p<0.05 were observed in all 13 genes with virologic response (Table 3), but none of these reached genome-wide significance. In order to further evaluate these nominally significant associations, gene-level burden testing analysis was conducted. None of the 13 burden tests was statistically significant ( Table 3).
All genetic effects on virologic response were stratified based on type of regimen to test associations specific to treatment with triple NRTI-based, NNRTI-based and PI-based regimens.
Associations between SNPs in these ADME genes and response to HAART were observed in almost all genes for the three regimen types. Again, no single association achieved genome-wide significance in any of the three regimen types. Gene-level associations were assessed for each type of treatment regimen to evaluate evidence of association after taking into account all of the variation in the region and accounting for correlation between the SNPs. No significant (P<0.05) burden tests were observed for the NNRTI-based regimen group (Table 4). However, three significant burden tests were observed in each of the PIbased regimen and the triple NRTI-alone groups.
CYP2B6, CYP2C19, and CYP2C9 were significantly associated with virologic response to HAART at the gene level among women on a PI-based regimen (p-values 0.049, 0.038 and 0.027, respectively). CYP2C9 was also associated with response to HAART in the triple NRTI-alone group (p=0.021) as were ADH1A and UGT1A1 (p-values 0.017 and 0.029, respectively). By chance, 2.5 significant burden tests would be expected compared to six observed.

Discussion
We evaluated the relationship between variation in 13 genes involved in the absorption, distribution, metabolism and elimination of ARVs and virologic response to an initial HAART regimen and observed twice as many gene-level burden tests with p-values less than 0.05 compared to what was expected. This suggests that variation in CYP2B6, CYP2C19, and CYP2C9 are associated with virologic response to PI-based regimens and variation in CYP2C9, ADH1A and UGT1A1 are associated with responses to regimens containing three NRTIs.
The lack of association overall (Table 3) is not surprising given the expected class-specific effects of these genes. CYP2B6 is only known to be involved directly in the metabolism of NNRTIs, though its activity is inhibited by the pharmacoenhancer, ritonavir [33]. The burden p-value for CYP2B6 was 0.049 for the PI based regimen and this association may likely be due to chance. CYP2C19 and CYP2C9 (burden p values of 0.038 and 0.027, respectively) are both known to be involved in the metabolism of PIs [34,35]. Genetic polymorphisms in both of these genes are known to markedly influence the clearance of other drugs [36].
ADH1A and UGT1A1 are the only genes known to be involved in the metabolism of NRTIs and the burden test p-values for these genes with triple NRTI regimens were 0.017 and 0.029, respectively. ADH1A and UGT1A1 are both involved in the metabolism of the NRTI abacavir [37]. SNPs in these regions could slow drug metabolism, possibly resulting in higher serum concentrations and more effective responses. CYP2C9 was also associated with NRTI-only treatment (burden test p value 0.021); as this gene was not previously known to be involved with NRTI metabolism, this could be a chance association. NRTI-only regimens are no longer recommended, but these results do provide insight into the metabolism of these drugs.
The major limitation of this study is the small sample size resulting in limited power to detect modest effects. Burden testing is a useful tool in this context to determine whether there is value in further follow-up of observed associations. Although some of the burden test results in our study may still be due to chance, the observation of six significant results suggests that one of these regions might be associated with virologic response to HAART and further study is warranted. The next steps in evaluating these associations are complicated by the fact that a region and not a SNP has been implicated. Thus for replication and extension of our findings, dense coverage of the regions would be needed and, optimally the population would also be female, as gender-specific effects cannot be ruled out. Fine mapping of the regions in the WIHS population that included rare variants is another option; if the true causal allele(s) are associated with a stronger effect, it might be possible to identify which of these six regions is truly associated with virologic response to HAART. In addition, some of these regions have poor coverage, particularly in the African American populations (Table 2), which likely limits the ability to detect significant associations for this group.
Though no genome wide assocations were detected, several gene level statistically significant associations were identified. Despite the low power of this study, at least one of these associations is likely to be true and not due to chance alone. Defining genetic factors that impact virologic response could lead to improvements in individualizing treatment for HIV-positive patients through genetic screening. Further studies are needed to more finely map these genetic regions harboring putative variants associated with virologic responses to specific HAART regimens.  Table 2 Description of the SNP coverage in the three racial/ethnic groups for the 13 genes studied.  Table 3 Gene level associations for response to HAART in all subjects, regardless of regimen type.  Table 4 Gene level associations between the 13 genes and response to HAART by type of regimen.