Pharmacogenetic Analysis of the Model-Based Pharmacokinetics of Five Anti-HIV Drugs: How Does This Influence the Effect of Aging?

Analysis of aging and pharmacogenetics (PGx) on antiretroviral pharmacokinetics (PKs) could inform precision dosing for older human HIV-infected patients. Seventy-four participants receiving either atazanavir/ritonavir (ATV/RTV) or efavirenz (EFV) with tenofovir/emtricitabine (TFV/FTC) provided PK and PGx information. Aging-PGx-PK association and interaction analyses were conducted using one-way analysis of variance (ANOVA), multiple linear regression, and Random Forest ensemble methods. Our analyses associated unbound ATV disposition with multidrug resistance protein (MRP)4, RTV with P-glycoprotein (P-gp), and EFV with cytochrome P450 (CYP)2B6 and MRP4 genetic variants. The clearance and cellular distribution of TFV were associated with P-gp, MRP2, and concentrative nucleoside transporters (CNTs), and FTC parameters were associated with organic cation transporters (OCTs) and MRP2 genetic variants. Notably, p16INK4a expression, a cellular aging marker, predicted EFV and FTC PK when genetic factors were adjusted. Both age and p16INK4a expression interacted with PGx on ATV and TFV disposition, implying potential dose adjustment based on aging may depend on genetic background.

Study Highlights

WHAT IS THE CURRENT KNOWLEDGE ON THE TOPIC?

✓ Previous studies determined an important role of aging on disposition of some ARVs. Gene polymorphisms, especially in metabolic enzymes, are related to variability in ARV PK, but their role in influencing the aging‐PK relationship is unknown.

WHAT QUESTION DID THIS STUDY ADDRESS?

✓ This study asked whether gene polymorphisms alter the interpretation of the effect of aging on (1) unbound ATV, RTV, and EFV clearance and (2) clearance and cellular distribution of TFV and FTC.

WHAT THIS STUDY ADDS TO OUR KNOWLEDGE

✓ This study updates the previous aging studies, shedding light on the significant effect of the cellular aging biomarker, p16INK4a, on unbound EFV disposition. We also present evidence of PGx‐aging interaction in ATV and TFV PKs, implying potential dose adjustment of these two drugs based on aging may be dependent on gene polymorphisms.

HOW THIS MIGHT CHANGE CLINICAL PHARMACOLOGY OR TRANSLATIONAL SCIENCE

✓ The p16INK4a might be an important biomarker to guide ATV dosing for older HIV‐infected patients. Evaluating the role of aging in dose adjustment of ATVs should take into account the presence of genetic variants in metabolic enzymes and drug transporters.

Gene polymorphisms play an important role in drug disposition and individualized therapy, as they can contribute to variability in pharmacokinetics (PKs), and may alter the effect of nongenetic factors. In HIV treatment, pharmacogenetic‐pharmacokinetic (PGx‐PK) relationships have been demonstrated for many antiretrovirals (ARVs) in the general population.1 Their roles in special populations, such as pregnancy2 and pediatrics,3 have also been evaluated. However, the PGx of ARVs in the older HIV‐infected population has not been well studied.

Aging is a trend in the HIV‐infected population.4 The estimated proportion of patients aged 50 years old, which is the general definition of “older” in the HIV field, was 42% in 2013,5 and is predicted to be 70% by the end of 2035.6 Aging could change drug absorption, distribution, metabolism, and elimination,7 which could alter efficacy and safety profiles. In our previous studies, we have investigated the effects of chronological age and aging surrogates on the PK of five selected ARVs.8, 9 Our analysis showed that the clearances of two renally eliminated drugs, tenofovir (TFV) and emtricitabine (FTC), demonstrated a significant association with creatinine clearance, which is strongly related to chronological age.9 In addition, we found that the cellular clearances of their active metabolites, tenofovir diphosphate (TFV‐DP) and emtricitabine triphosphate (FTC‐TP), are associated with the expression of a cellular aging biomarker, p16INK4a.9 On the other hand, three hepatically metabolized drugs, atazanavir (ATV), ritonavir (RTV), and efavirenz (EFV), did not show statistically significant relationships with aging, in either total or unbound drug disposition,8 although the sample size for these drugs was smaller than that of TFV/FTC.

Identification of influential single‐nucleotide polymorphisms (SNPs) could provide PGx information for ARV dosing adjustment, and improve statistical analyses of aging by reducing variability in PK. Additionally, gene polymorphisms can modify the effect of aging on PK, and evidence for such interactions has been shown for non‐ARV drugs. For example, less profound changes in systemic clearance of omeprazole between younger and older patients was observed in cytochrome P450 (CYP)2C19 poor metabolizers than in extensive metabolizers.10 This implies that dose adjustment in older patients might be less important given the dysfunctional CYP2C19 genotype. In another example, an eightfold increase in the mean dose‐adjusted plasma concentration of venlafaxine was seen in older compared with younger patients among the CYP2D6 poor metabolizers, but not in the other genotype groups.10 This suggests that dose adjustment might be necessary in the presence of CYP2D6 dysfunction. Therefore, investigating aging in the context of PGx or aging‐PGx interaction is essential to provide precise treatment for older patients.

In addition to the aging‐PK associations, as outlined in our previous publications, our studied drugs also serve as reasonable models for aging‐PGx‐PK investigation, because their disposition profiles cover a broad spectrum of metabolic enzymes and drug transporters. ATV and RTV are CYP3A substrates,11, 12 and served as two model drugs to probe these enzymes. In addition, the upstream CYP3A expression regulator, pregnane X receptor (PXR), affects ATV clearance.13, 14 Uridine 5′‐diphosphate glucuronosyltransferase (UGT)1A is responsible for ATV glucuronidation and is inhibited by ATV.11 EFV has a complex metabolic profile, which involves CYP2B, CYP2A, CYP3A, CYP1A, and UGT.15 It is also an inducer or inhibitor of CYP2B, CYP3A, and CYP2C.16 Transporters in the gut, such as P‐glycoprotein (P‐gp) and multidrug resistance proteins (MRPs), are presumably important globally to ARVs, which are mostly orally administrated.1 Organic anion‐transporting polypeptides (OATPs), which are enriched in the liver, could be relevant to ATV, RTV, and EFV, and OATs or organic cation transporters (OCTs) in the kidney are important to TFV or FTC.17 Nucleoside transporters, concentrative nucleoside transporters (CNTs) and equilibrative nucleoside transporters, are believed to transport nucleoside analogs,1, 17 which could be crucial to the disposition of TFV and FTC and the cellular accumulation of their active metabolites.

In the current study, we aim to re‐evaluate the effect of aging on PK on five studied drugs in the context of PGx and investigate the aging‐PGx interaction in PK. Using a candidate gene approach to select important SNPs, we also aim to provide new hypotheses for further PGx investigations of ARVs. For phenotype, we have focused on the clearance of TFV, FTC and unbound ATV, RTV and EFV, as well as cellular clearance of TFV‐DP and FTC‐TP, using the associated individual parameters derived from the previous established population PK models.

MATERIALS AND METHODS

Clinical trial

Study conduct has been described in previous publications.8, 9 Briefly, HIV‐infected patients were enrolled from the University of North Carolina HealthCare Infectious Diseases Clinic (Chapel Hill, NC) and the Cone Health Regional Center for Infectious Diseases (Greensboro, NC). All the subjects received steady‐state tenofovir disoproxil fumarate (TDG)/FTC 300/200 mg with either EFV 600 mg or ATV/RTV 300/200 mg daily. Blood samples were collected either sparsely at four time points (predose, and then 2 h, 4–6 h, and 10–14 h postdose) or intensively at 11 time points for PK analysis. An extra blood sample was collected for those who consented to be involved into the PGx analysis. The study protocol was approved by the institutional review boards of both institutions (Clinicaltrials.gov NCT01180075).

Genotyping

Samples for this test were obtained from patients who consented to genotyping assays at screening. Genomic DNA was extracted from whole blood using QIAamp DNA Blood Mini Kit (QIAGEN). Genotyping was conducted for selected metabolic enzymes and drug transporters using DMET Plus microarray (Affymetrix, CA). Samples with a call rate <95% were excluded from analysis. The SNPs were excluded if its minor allele frequency was <4% and its missing genotypes were >30%. Deviations from expected proportions based on Hardy‐Weinberg equilibrium (exact test, α > 0.05) and linkage disequilibrium (R2 > 0.2) analysis were performed in PLINK software version 1.07.18

Phenotypes

The phenotypes in this study were the post hoc PK parameters. The PK analysis was performed as previously described.8, 9 Briefly, drug concentrations were measured in the University of North Carolina Center for AIDS Research Clinical Pharmacology and Analytical Chemistry Laboratory. Total and intracellular drug concentrations were measured by validated liquid‐chromatography tandem mass spectrometry. Unbound EFV, ATV, and RTV concentrations were analyzed using rapid equilibrium dialysis followed by liquid‐chromatography tandem mass spectrometry. Population PK models were developed for each drug using nonlinear mixed effects modeling program NONMEM version 7.3 (ICON Development Solutions, Hanover, MD). For the current study, individual PK parameters of interest included clearance of the unbound drug for EFV, ATV, and RTV, clearance of the parent drug for TFV and FTC, and the rate constant of intracellular metabolite conversion and elimination for TFV‐DP and FTC‐TP. These estimates were obtained using the Bayesian post hoc estimation method based on the final models, including the covariate effects of creatinine clearance on TFV and FTC clearance and p16INK4a expression on TFV‐DP and FTC‐TP disposition.

Statistical analysis

A candidate gene approach was used to perform this analysis. The candidate genes were chosen if shown to be relevant to studied drug disposition pathways. Variables were log‐transformed to achieve normality when appropriate. In the univariate analysis, the association analysis between genotypes and PK parameters was conducted using one‐way analysis of variance (ANOVA). Then, SNPs with a nominally significant P value were included in multivariate analysis using the stepwise selection procedure in both forward and backward directions. Model selection was based on Akaike information criteria, involving demographic covariates of interest, including chronological age, p16INK4a expression level, race, sex, and body mass index (BMI). Missing genotypes were imputed using IMPUTE2 with the 1000 Genomes Phase III panel.19, 20 For race, patients who were other, unknown, and nonreported races were grouped as one category (“others”) due to low counts. Additionally, the effect of creatinine clearance and treatment arms were evaluated for TFV and FTC. Comparison of demographics between two treatment arms was conducted using Mann‐Whitney U test or Kruskal‐Wallis test. These statistical analyses were performed using R version 3.2 (R‐project.org). A P value < 0.05 was considered statistically significant.

To investigate the interaction between aging and gene polymorphism, we used distinct linear regression models (R function lm, package stats version 3.4.0; www.r-project.org)21 to determine if there is an effect on drug clearance due to the interaction between chronological age and gene polymorphisms (SNPs), and between p16INK4a expression and gene polymorphisms. An F‐test (R function ANOVA, package stats version 3.4.0)21 was used to compare the model without interaction effects to the one with interaction effects. Only significant predictors from the best fitting model were retained in the final model, with the cutoff for significance being P < 0.05 for each predictor. The factors involved in each interaction effect were always included in the model. Multiple testing corrections were performed across all drugs using the Benjamini‐Hochberg method with a false discovery rate of q < 0.1 (R function p.adjust, package stats version 3.4.0).21

In order to determine variable interactions and variable importance of all available variables, we used the Random Forest ensemble method (R function randomForest(ntree = 5000), package randomForest version 4.6‐12)22 for each drug, with the PK parameter as the response, after removing all incomplete cases. We used the percent increase in mean squared error to determine the three most important predictors for each drug.

The effects of significant SNPs on corresponding proteins were evaluated based on literature and database (e.g., dbSNP, PharmGKB) searching. The SIFT and PolyPhen algorithms were used to make predictions for SNPs with no current information on how they influence the protein. If still no information is available, those SNPs in linkage disequilibrium with r2 >0.8 were also investigated. Protein function prediction was conducted using SNPnexus.23, 24, 25 Linked SNPs were obtained through HaploReg version 4.1.26

RESULTS

Population

In 91 patients enrolled in the original clinical study, 74 consented to the genetic tests; 25 were in the ATV/RTV/TFC/FTC arm and 49 in the EFV/TFV/FTC arm. The median age of the 74 patients was 48 years. The self‐reported race includes the following categories: white (24/74); Black/African American (45/74); other (3/74); unknown (1/74), and not reported (1/74). Fifty‐two percent of the patients were men. Three patients demonstrated the Fried frailty phenotype.27 The characteristics between the two treatment arms were not significantly different. The studied PK parameters included apparent unbound drug clearance (CLu/F) of ATV, RTV, and EFV, apparent clearance (CL/F) of TFV and FTC, and clearance of TFV and FTC metabolites (TFV‐DP and FTC‐TP, respectively) into (CLin) and out of (CLout) of peripheral blood mononuclear cells. In total, 22 genes were selected for investigation. Patient demographics, PK parameters, and studied genes and SNPs are summarized in Table 1.

Table 1

Overview of demographics of studied population and the pharmacokinetics and pharmacogenetics of the studied drugs

		ATV	RTV	EFV	TFV	FTC
No. of subjects		25		49	74
Age (years)		49 (24, 61)		48 (22, 73)	48 (22, 73)
Race	African	15 (60%)		30 (61%)	45 (61%)
	American
	White	8 (32%)		16 (33%)	24 (32%)
	Others	2 (8%)		3 (6%)	5 (7%)
Sex	Female	10 (40%)		12 (24%)	22 (30%)
	Male	15 (60%)		37 (76%)	52 (70%)
CrCL (L/h)		96.6 (66.8, 228)		112 (43.3, 200)	109 (43.3, 228)
BMI (kg/m2)		29.8 (20.2, 40.4)		27.1 (17.6, 44.3)	28.6 (17.6, 44.3)
Frailty status	Nonfrail	16 (64%)		37 (76%)	53 (72%)
	Pre‐frail	8 (32%)		10 (20%)	18 (24%)
	Frail	1 (4%)		2 (4%)	3 (4%)
p16INK4a expression a		2.20 (0.158, 3.91)		1.92 (‐1.14, 2.81)	2.01 (‐1.14, 3.91)

ATV, atazanavir; BMI, body mass index; CL/F, apparent clearance; CLin, clearance into the cell; CLout, clearance out of the cell; CLu/F, apparent unbound drug clearance; CrCL, creatinine clearance; CYP, cytochrome P450; EFV, efavirenz; F, undetermined bioavailability; FTC, emtricitabine; FTC‐TP, emtricitabine triphosphate; PK, pharmacokinetic; PXR, pregnane X receptor; RTV, ritonavir; SNP, single‐nucleotide polymorphism; TFV, tenofovir; TFV‐DP, tenofovir diphosphate; UGT, uridine 5′‐diphosphate glucuronosyltransferase.

Clearance is in liter/hour (L/h).

aThe expression of p16INK4a is presented as logarithm (base 2) of gene expression. Continuous variables are in median (range). Categorical variables are in count (percentage).

Linear regression analysis

The P values for univariate analysis are shown in Table 2. Twenty‐one SNPs were found to be significantly associated with the PK parameters. The SNP basic information, genotype frequency, effects on phenotypes, and SNP function are summarized in Table 3.15, 28, 29, 30, 31 No significant SNPs were found for the TFV‐DP CLout. The FTC‐TP CLin showed the same statistical results as the FTC CL/F, because during the model construction the FTC‐TP CLin was described as a fixed fraction of FTC CL/F.9 Therefore, TFV‐DP CLout and FTC‐TP CLin were excluded from further analysis. Associations between demographics and PK parameters were consistent with the previous population PK analysis. Briefly, BMI was associated with RTV clearance. Age and creatinine clearance were correlated, and shown to be related to CL/F of both TFV and FTC. The p16INK4a expression was associated with FTC CL/F and FTC CLout.

Table 2

The P values in the univariate analysis

Factors	ATV CLu/F	RTV CLu/F	EFV CLu/F
Age	0.74	0.64	0.55
p16 expression	0.20	0.43	0.57
Race	0.40	0.95	0.22
BMI	0.66	0.025*	0.64
Sex	0.10	0.75	0.07
Frailty status	0.074	0.21	0.69
CrCL	–	–	–
Treatment arm	–	–	–
SNP	SLCO1B1_c.571T>C: 0.036 ABCC4_c.*879T>C: 0.040	ABCB1_c.3435C>T: 0.0021 CYP1A2 * 1C_‐3860G>A: 0.011 CYP2B6_14593C>G: 0.047	CYP2B6 * 6_15631G>T: 0.00025 ABCC4_c.3348A>G/Del: 0.0056 ABCC4_c. *311G>A: 0.0055 CYP2B6 * 4_18053A>G(K262R): 0.013 ABCC4_c. *38T>G: 0.013

ATV, atazanavir; BMI, body mass index; CL/F, apparent clearance; CLin, clearance into the cell; CLout, clearance out of the cell; CLu/F, apparent unbound drug clearance; CrCL, creatinine clearance; EFV, efavirenz; F, undetermined bioavailability; FTC, emtricitabine; FTC‐TP, emtricitabine triphosphate; RTV, ritonavir; SNP, single‐nucleotide polymorphism; TFV, tenofovir; TFV‐DP, tenofovir diphosphate.

* P < 0.05; *** P < 0.001.

Table 3

Summary of the SNPs showing statistical significance in the univariate analysis

Studied drug	dbSNP RS	Associated gene/protein	c.	Base change	Function	Genotype frequency	Phenotype	Previously reported phenotype	Physiological reasonable?	Physiological reason
ATV	rs1059751	ABCC4/MRP4	13	879T>C	3'‐UTR	T/T+T/C: 23 (92%) C/C: 2 (8%)	↑ CL/F	Decreased function	Y	Decreased efflux from liver to system
	rs4149057	SLCO1B1/OATP1B1	12	571T>C	L191L	T/T+T/C: 21 (84%) C/C: 4 (16%)	↑ CL/F	–	–	–
RTV	rs1045642	ABCB1/P‐gp	7	3435C>T	I1145I	C/C: 14 (56%) C/T+T/T: 11 (44%)	↓ CL/F	Decreased expression	Y	Increased bioavailability
	rs2069514	CYP1A2/CYP1A2	15	‐3860G>A	Promoter	G/G+G/A: 22 (88%) A/A: 3 (12%)	↓ CL/F	1C,a decreased function	Y	Decreased metabolism
	rs4803418	CYP2B6/CYP2B6	19	14593C>G	Intronic	C/C: 19 (76%) C/G: 6 (24%)	↑ CL/F	1C,a increased function	Y	Increased metabolism
EFV	rs3745274b	CYP2B6/CYP2B6	19	15631G>T	Q172H	G/G+G/T: 41 (84%) T/T: 6 (12%) Missing: 2 (4%)	↓ CL/F	6,a decreased function	Y	Decreased metabolism
	rs2279343	CYP2B6/CYP2B6	19	18053A>G	K262R	A/A+A/G: 39 (80%) G/G: 7 (14%) Missing: 3 (6%)	↓ CL/F	4,a decreased function	Y	Decreased metabolism
	rs1751034	ABCC4/MRP4	13	3348A>G/Del	K1116K/X	A/A+A/G: 45 (92%) G/G: 4 (8%)	↓ CL/F	Decreased function	Y	Increased bioavailability
	rs4148551	ABCC4/MRP4	13	311G>A	3’‐UTR	G/G: 11 (22%) G/A+A/A: 38 (78%)	↑ CL/F	–	–	–
	rs3742106	ABCC4/MRP4	13	38T>G	3’‐UTR	T/T+T/G: 46 (94%) G/G: 3 (6%)	↓ CL/F	Decreased function	Y	Increased bioavailability
TFV	rs3213619	ABCB1/P‐gp	7	‐129T>C	5’‐UTR	T/T: 61 (82%) T/C: 13 (18%)	↓ CL/F	Decreased function	Y	Increased bioavailability
	rs2242046	SLC28A1/CNT1	15	1561G>A	D521N	G/G: 49 (66%) G/A+A/A: 19 (26%) Missing: 6 (8%)	↓ CL/F	Increased function	Y	Decreased renal excretion
	rs8023604	SLC28A2/CNT2	15	531T>C	F177F	T/T+T/C: 58 (79%) C/C: 15 (20%) Missing: 1 (1%)	↑ CL/F	–	–	–
TFV‐DP	rs2235013	ABCB1/P‐gp	7	1725+38G>A	Intronic	G/G: 24 (27%) G/A+A/A: 50 (73%)	↓ CLin	–	–	–
	rs10519020	SLC28A2/CNT2	15	734G>C	S245T	G/G+G/C: 72 (97%) C/C: 2 (3%)	↑ CLin	–	Na	Presumably decreased cellular uptake
	rs17222723	ABCC2/MRP2	10	3563T>A	V1188E	T/T: 64 (86%) T/A: 10 (14%)	↑ CLin	Decreased functionc	Y	Decreased efflux
FTC	rs668871	SLC22A3/OCT3	6	360C>T/G	R120R	C/C: 14 (19%) C/T+T/T: 57 (77%) Missing: 3 (4%)	↑ CL/F	Increased expressiond	Y	Increased renal excretion
	rs316003	SLC22A2/OCT2	6	1506A>G	V502V	A/A+A/G: 60 (81%) G/G: 14 (19%) Missing: 1 (1%)	↑ CL/F	–	–	–
	rs8187707	ABCC2/MRP2	10	4488C>T	H1496H	C/C: 67 (91%) C/T: 7 (8%)	↑ CL/F	Decreased functionc	N	Presumably decreased renal excretion or increased bioavailability
FTC‐TP	rs2273697	ABCC2/MRP2	10	1249G>A	V417I	G/G: 45 (61%) G/A+A/A: 29 (39%)	↑ CLout/F	Increased function	Y	Increased efflux
	rs2305367	SLC28A1/CNT1	15	1149A>G	K383K	A/A+A/G: 44 (60%) G/G: 27 (36%) Missing: 3 (4%)	↑ CLout/F	–	–	–

ATV, atazanavir; c., chromosome number; CL/F, apparent clearance; CLin, clearance into the cell; CLout, clearance out of the cell; CNT, concentrative nucleoside transporter; CYP, cytochrome P450; dbSNP, database‐single‐nucleotide polymorphism; EFV, efavirenz; F, undetermined bioavailability; FTC, emtricitabine; FTC‐TP, emtricitabine triphosphate; MRP, multidrug resistance protein; OATP, organic anion‐transporting polypeptide; P‐gp, P‐glycoprotein; RS, reference single‐nucleotide polymorphism cluster ID; SNP, single‐nucleotide polymorphism; TFV, tenofovir; TFV‐DP, tenofovir diphosphate; UTR, untranslated region.

The finding is inconsistent with the predicted phenotypic effects.

The rs8187710 (r2 = 0.98 with rs17222723 and rs8187707) was shown to decrease MRP2 function (ref. 30 and 31).

The rs555754 (r2 = 0.99 with rs668871) was shown to increase OCT3 messenger RNA expression (ref. 29).

The rs3745274 is a well‐identified SNP influencing EFV elimination (ref. 15 and 28).

ABCB1_c.‐129T>C: 0.013

SLC28A1_c.1561G>A: 0.015

SLC28A2_c.531T>C: 0.0060

ABCB1_c.1725+38A>G: 0.0064

SLC28A2_c.734G>C: 0.0090

ABCC2_c.3563T>A: 0.038

SLC22A3_c.360C>T>G: 0.032

SLC22A2_c.1506G>A: 0.018

ABCC2_c.4488C>T: 0.020

ABCC2_c.1249G>A: 0.0046

SLC28A1_c.1149A>G: 0.038

In the multiple linear regression, ABCC4_c.*879T>C (rs1059751) was associated with increased ATV CLu/F, with no demographics being statistically significant. Decreased RTV CLu/F was associated with the mutation at ABCB1_c.3435 (rs1045642). BMI remained significant with RTV CLu/F. In addition, race was found to be predictive of RTV CLu/F; African Americans had significantly lower apparent clearance of RTV than the other populations. Decreased EFV CLu/F was related to variants CYP2B6*6_15631G>T (rs3745274), ABCC4_c.3348A>G/Del (rs1751034), and ABCC4_c.*38T>G (rs3742106). Notably, after adjusting for the effects of these SNPs, p16INK4a expression was found to be significantly associated with EFV CLu/F. A doubling of p16INK4a expression was shown to be related to 315 L/h (31% of the median) decrease in EFV unbound CLu/F.

The polymorphisms at ABCB1_c.‐129 (rs3213619) and SLC28A1_c.1561 (rs2242046) were significantly associated with decreased TFV CL/F, whereas SLC28A2_c.531T>C (rs8023604) showed the opposite effect. Creatinine clearance remained significant, whereas age was removed from the model. Sex was related to TFV CL/F after adjusting for gene polymorphism, with women demonstrating lower TFV clearance than men. The ABCB1_c.1725+38G>A (rs2235013), SLC28A2_c.734G>C (rs10519020), and ABCC2_c.3563T>A (rs17222723) were retained in the model for TFV‐DP CLin. For FTC CL/F, the significant predictors were SLC22A3_c.360C>T>G (rs668871), SLC22A2_c.1506A>G (rs316003), ABCC2_c.4488C>T (rs8187707), creatinine clearance, and p16INK4a expression; the effect of age was eliminated. The ABCC2_c.1249G>A (rs2273697) and p16INK4a expression were the two significant predictors in the FTC CLout linear model. Results from the multiple linear regressions are summarized in Table 4.

Table 4

Stepwise linear regression analysis

PK parameters	Factor	P value	Effect size
ATV CLu/F	ABCC4_c.*879T>C (MRP4)	0.040*	151.4
RTV CLu/F	ABCB1_c.3435C>T (P‐gp)	<0.001***	−855.6
	BMI	<0.001***	70.21
	Caucasiana	0.0027**	603.1
	Othersa	0.0036**	756.9
EFV CLu/F	CYP2B6*6_15631G>T (CYP2B6)	<0.001***	−710.2
	ABCC4_c.3348A>G/Del (MRP4)	0.0020**	−768.7
	ABCC4_c.*38T>G (MRP4)	<0.001***	−845.8
	p16INK4 a expression	0.0083**	−315.2
TFV CL/F	ABCB1_c.‐129T>C (P‐gp)	<0.001***	−15.22
	SLC28A1_c.1561G>A (CNT1)	0.0047**	−12.48
	SLC28A2_c.531T>C (CNT2)	<0.001***	19.18
	CrCL	<0.001***	0.3220
	Sex – female	0.049*	−6.660s
TFV‐DP CLin	ABCB1_c.1725+38G>A (P‐gp)	0.0031**	−0.0006431
	SLC28A2_c.734G>C (CNT2)	0.028*	0.0026217
	ABCC2_c.3563T>A (MRP2)	0.033*	0.0006769
FTC CL/F	SLC22A3_c.360C>T>G (OCT3)	0.013*	3.882
	SLC22A2_c.1506A>G (OCT2)	0.019*	3.860
	ABCC2_c.4488C>T (MRP2)	0.0063**	6.026
	CrCL	<0.001***	0.07060
	p16INK4 a expression	0.034**	−2.006
FTC‐TP CLout	ABCC2_c.1249G>A (MRP2)	0.022*	0.01012
	p16INK4 a expression	0.0090**	0.008336

ATV, atazanavir; BMI, body mass index; c., chromosome number; CL/F, apparent clearance; CLin, clearance into the cell; CLout, clearance out of the cell; CLu/F, apparent unbound drug clearance; CNT, concentrative nucleoside transporter; CrCL, creatinine clearance; EFV, efavirenz; F, undetermined bioavailability; FTC, emtricitabine; FTC‐TP, emtricitabine triphosphate; MRP, multidrug resistance protein; OCT, organic cation transporter; P‐gp, P‐glycoprotein; PK, pharmacokinetic; RTV, ritonavir; TFV, tenofovir; TFV‐DP, tenofovir diphosphate.

*P < 0.05; **P < 0.01; ***P < 0.001.

Reference race: African American; “others” includes other races, unknown races, and unreported races.

Note: Only factors that are retained in the final linear regression model are shown in this table. Drug clearance was modeled on the log scale. The expression of p16INK4a is presented as logarithm (base 2) of gene expression. Protein names related to the genes are presented in the parentheses.

Interaction analysis

To analyze interactions, we used two distinct linear regression models to determine if there is an interaction effect between chronological age and genetic polymorphisms, and between p16INK4a expression and genetic polymorphisms. Most of the significant SNPs in the previous multiple linear regressions were confirmed in this analysis, which was without other demographic variables (Table 5). Notably, we found significant evidence of the interaction effects of age and ABCC4_c.*879T>C (rs1059751) with ATV CLu/F. This effect was also found in the model with p16INK4a expression. In addition, both models showed very strong significant evidence of the interaction between age or p16INK4a expression and SLC28A2_c.734G>C (rs10519020) with TFV‐DP CLin.

Table 5

Aging‐PGx interaction‐focused linear regression analysis

PK parameters	Interaction model with	Factor	P value	Q value	Effect size
ATV CLu/F	Age	ABCC4_c.*879T>C	3.60E‐02	4.02E‐02	−1.33E+03
		Age:ABCC4_c.*879T>C	2.07E‐02	2.77E‐02	3.29E+01
	p16INK4a expression	p16INK4a expression:ABCC4_c.*879T>C	1.44E‐02	2.28E‐02	1.18E+02
RTV CLu/F	Age	ABCB1_c.3435C>T	3.47E‐03	7.32E‐03	−6.29E+02
		CYP2B6*1C_14593C>G	5.75E‐03	1.03E‐02	6.83E+02
	p16INK4a expression	ABCB1_c.3435C>T	1.63E‐03	6.18E‐03	−6.99E+02
		CYP2B6*1C_14593C>G	2.07E‐02	2.62E‐02	5.71E+02
EFV CLu/F	Age	CYP2B6*6_15631G>T	3.77E‐04	1.19E‐03	−6.92E‐01
		ABCC4_c.*38T>G	4.77E‐04	1.29E‐03	−8.73E‐01
		ABCC4_c.3348A>G/Del	2.58E‐02	3.06E‐02	−4.94E‐01
	p16INK4a expression	p16INK4a expression	6.60E‐03	1.27E‐02	−2.41E‐01
		CYP2B6*6_15631G>T	3.04E‐04	1.44E‐03	−6.86E‐01
		ABCC4_c.*38T>G	2.48E‐04	1.44E‐03	−8.66E‐01
		ABCC4_c.3348A>G/Del	2.79E‐03	6.74E‐03	−8.16E‐01
TFV CL/F	Age	Age	2.55E‐07	1.62E‐06	−1.11E+00
		SLC28A2_c.531T>C	5.95E‐03	1.03E‐02	1.52E+01
	p16INK4a expression	SLC28A2_c.531T>C	3.94E‐02	4.68E‐02	1.34E+01
		SLC28A1_c.1561G>A	1.63E‐02	2.38E‐02	−1.47E+01
TFV‐DP CLin	Age	SLC28A2_c.734G>C	8.28E‐10	1.57E‐08	3.58E‐02
		Age: SLC28A2_c.734G>C	8.45E‐09	8.03E‐08	−6.20E‐04
		ABCB1_c.1725+38G>A	2.19E‐02	2.77E‐02	−3.89E‐04
	p16INK4a expression	SLC28A2_c.734G>C	3.29E‐11	6.26E‐10	1.68E‐02
		p16INK4a expression: SLC28A2_c.734G>C	4.99E‐09	4.74E‐08	−6.14E‐03
FTC CL/F	Age	Age	2.08E‐04	7.91E‐04	−1.65E‐01
		ABCC2_c.4488C>T	1.45E‐04	6.90E‐04	6.29E+00
		SLC22A2_c.1506A>G	2.45E‐03	5.83E‐03	3.72E+00
		SLC22A3_c.360C>T/G	1.07E‐02	1.70E‐02	3.11E+00
	p16INK4a expression	p16INK4a expression	1.86E‐02	2.52E‐02	−1.81E+00
		ABCC2_c.4488C>T	2.17E‐03	6.74E‐03	6.32E+00
		SLC22A2_c.1506A>G	2.84E‐03	6.74E‐03	4.00E+00
		SLC22A3_c.360C>T/G	9.33E‐03	1.61E‐02	3.65E+00
FTC‐TP CLout	Age	ABCC2_c.1249G>A	1.63E‐02	2.38E‐02	1.13E‐02
	p16INK4a expression	p16INK4a expression	6.69E‐03	1.27E‐02	9.23E‐03

ATV, atazanavir; c., chromosome number; CL/F, apparent clearance; CLin, clearance into the cell; CLout, clearance out of the cell; CLu/F, apparent unbound drug clearance; EFV, efavirenz; F, undetermined bioavailability; FTC, emtricitabine; FTC‐TP, emtricitabine triphosphate; PGx, pharmacogenetics; RTV, ritonavir; TFV, tenofovir; TFV‐DP, tenofovir diphosphate.

The P values, Q values, and effect sizes shown from interaction models containing gene polymorphisms (single‐nucleotide polymorphisms [SNPs]) and age as factors or SNPs and p16INK4a expression as factors. Only factors that were significant after a multiple testing correction at the false discovery rate of q < 0.1 are shown. Drug clearance was modeled on the log scale for drug EFV. The expression of p16INK4a is presented as logarithm (base 2) of gene expression.

A colon between two variables indicates interaction between the two.

In addition to the analysis with traditional statistical approaches, we performed machine learning analysis to evaluate the variable interactions and the overall importance of each of the variables with the PK parameters. Although the Random Forest method did not find consistent interactions, we used it to determine variable importance. The three most important predictors for ATV CLu/F were ABCC4_c.*879T>C (rs1059751), SLCO1B1_c.571T>C, (rs4149057), and sex of the individual. For RTV CLu/F, they were ABCB1_c.3435C>T (rs1045642), CYP1A2*1C_‐3860G>A (rs2069514), and BMI. For EFV CLu/F, the three most important variables were CYP2B6*6_15631G>T (rs3745274), ABCC4_c.3348A>G/Del (rs1751034), and ABCC4_c. 311G>A (rs4148551). For TFV they were creatinine clearance, age, and SLC28A1_c.1561G>A (rs2242046). For TFV‐DP CLin, they were ABCB1_c.1725+38G>A (rs2235013), sex, and BMI. For FTC CL/F, they were creatinine clearance, p16INK4a expression, and age. For FTC‐TP CLout, the three most important predictors were p16INK4a expression, ABCC2_c.1249G>A (rs2273697), and SLC28A1_c.1149A>G (rs2305367). The three most important predictors for each drug are listed in Table 6.

Table 6

Top three predictors for drug PK in the Random Forest analysis

PK parameters	Top three predictors
ATV CLu/F	ABCC4_c.*879T>C (MRP4)
	SLCO1B1_c.571T>C (OATP1B1)
	Sex
RTV CLu/F	ABCB1_c.3435C>T (P‐gp)
	CYP1A2*1C_‐3860G>A (CYP1A2)
	BMI
EFV CLu/F	CYP2B6*6_15631G>T (CYP2B6)
	ABCC4_c.3348A>G/Del (MRP4)
	ABCC4_c. 311G>A (MRP4)
TFV CL/F	SLC28A1_c.1561G>A (CNT1)
	Age
	CrCL
TFV‐DP CLin	ABCB1_c.1725+38G>A (P‐gp)
	Sex
	BMI
FTC CL/F	Age
	CrCL
	p16INK4a expression
FTC‐TP CLout	ABCC2_c.1249G>A (MRP2)
	SLC28A1_c.1149A>G (CNT1)
	p16INK4a expression

ATV, atazanavir; BMI, body mass index; c., chromosome number; CL/F, apparent clearance; CLin, clearance into the cell; CLout, clearance out of the cell; CLu/F, apparent unbound drug clearance; CNT, concentrative nucleoside transporter; CrCL, creatinine clearance; CYP, cytochrome P450; EFV, efavirenz; F, undetermined bioavailability; FTC, emtricitabine; FTC‐TP, emtricitabine triphosphate; MRP, multidrug resistance protein; OATP, organic anion‐transporting polypeptide; P‐gp, P‐glycoprotein; PK, pharmacokinetic; RTV, ritonavir; TFV, tenofovir; TFV‐DP, tenofovir diphosphate.

Drug clearance was modeled on the log scale for drug EFV. The expression of p16INK4a is presented as logarithm (base 2) of gene expression. Protein names related to the genes are presented in the parentheses.

DISCUSSION

In this study, we investigated the effects of aging and gene polymorphisms on key PK parameters of five ARVs. Some participant characteristics were significant predictors after adjusting for genetic factors, demonstrating the importance and potential of accounting for PGx in PK analysis. The results from the interaction analysis brought to light interactions between genotype and aging (both chronological and cellular p16INK4a expression). The multiple analyses methods used reinforced several of our findings.

Our analysis showed that increased apparent clearance of unbound ATV is associated with ABCC4_c.*879T>C. ATV is highly metabolized in the liver; this variant might lead to decreased efflux32 of ATV out of the hepatic cells, thus increased ATV metabolism due to accumulation. There is a nominal association with rs4149057 in SLCO1B1, the gene that also showed marginal significance in a genomewide association study of ATV PK.33 Possibly due to limited sample size, our study was not able to detect the UGT1A1 genotype, such as UGT1A1 *28, as a significant predictor of ATV PK.

The PGx profile of RTV is poorly studied. In our analysis, we found that ABCB1_c.3435C>T is significantly associated with decreased CLu/F, as previously reported.34 RTV is a substrate of P‐pg;12 the presence of variant might increase the bioavailability of RTV, thus decreases the CLu/F. In the univariate analysis, CYP1A2*1C_‐3860G>A (CYP1A2 *1C, decreased protein function35) is associated with decreased RTV CLu/F, whereas CYP2B6_14593C>G (CYP2B6 *1C, conflicting protein function28, 35) shows the opposite effect. BMI remains a significant predictor and is associated with higher RTV CLu/F. Additionally, we found that race plays a role, with African Americans showing the lower RTV clearance than other races. However, the effect of “other” races is not interpretable due to the lack of information of the actual racial background. We found no association between race and the three SNPs mentioned above (data not shown), implying that race might predict RTV CLu/F independent of these three SNPs. Our finding suggests that RTV drug disposition is influenced by several factors in a mixed manner, potentially explaining the RTV observed high‐PK variability.

Our analyses replicated the demonstrated effect of CYP2B6*6_15631G>T (rs3745274) to decrease EFV clearance. In addition, we found ABCC4_c.3348A>G and ABCC4_c.*311G>A, both of which decrease MRP4 function,37, 38, 39 were associated with decreased EFV CLu/F. This might be due to increased bioavailability from MRP4 dysfunction. Interestingly, we found that after accounting for these SNPs, decreased EFV CLu/F was significantly related to higher p16INK4a expression, suggesting that elimination of EFV is associated with cellular aging. The underlying mechanism still needs to be elucidated, although previous work demonstrated that the expression of CYP2B6 is negatively regulated through a CDK2 signaling pathway involving p16INK4a in a hepatic carcinoma cell line, HepG2.40 Our analysis shows that EFV CLu/F decreases by 31% of the median with twofold increase in p16INK4a expression. Further analysis with larger sample size is needed to validate this finding and inform the potential application of this biomarker to dosing guidance in older HIV‐infected patients.

The significant SNPs in the univariate analysis independently predict TFV and FTC CL/F, given that all of them are retained in the final model. We found that TFV CL/F is associated with three variants on P‐gp, CNT1, and CNT2, respectively. The ABCB1_c.‐129T>C (rs3213619) is shown to decrease CL/F, which might be associated with increased bioavailability through the gut due to decrease P‐gp function.41 The SLC28A1_c.1561G>A (rs2242046) also decreases CL/F, and this may be related to the decreased renal excretion of TFV resulting in increased CNT1 function,42 which transports more drugs back to the nephrons. There is no information about SLC28A2_c.531T>C (rs8023604) function in the literature, but we hypothesize this SNP might decrease CNT2 function, in turn increasing drug excretion to the urine. The three SNPs associated with FTC elimination all increase CL/F. The SLC22A3_c.360C>T might be associated with increased OCT3 expression,29 which increases renal excretion of drugs. The ABCC2_c.4488C>T was implicated in one previous study with decreased MRP2 function,30, 31 presumably resulting in decreased renal excretion or increased bioavailability, whereas our results indicate the opposite. Additional work will be needed to clarify this finding; no functional information for the SLC22A2_c.1506A>G polymorphism is available. In terms of demographics, creatinine clearance is a significant predictor for both TFV and FTC CL/F, and, expectedly, the effect of age was eliminated in multivariate analysis, due to the collinearity of the two variables. However, p16INK4a expression was significantly associated with FTC CL/F, but not TFV CL/F; potentially, pathways sensitive to cellular aging might play a role in FTC elimination.

Cellular uptake of TFV‐DP (CLin) is significantly associated with three SNPs on P‐gp, CNT2, and MRP2 genes, respectively. TFV is not a substrate of P‐gp, thus this association might be due to the bioavailability of tenofovir disoproxil fumarate, the prodrug of TFV, during the absorption phase rather than the cellular distribution. The SLC28A2_c.734G>C (rs10519020) is shown to increase TFV‐DP CLin, implying this SNP may increase CNT2 function. However, there is no known information about this SNP, and the function prediction is the opposite of our finding. The ABCC2_c.3563T>A (rs17222723) is related to increased CLin, which makes some sense given that this SNP is associated with decreased function of MRP2.30, 31 There are some discordant findings regarding the relationship between TFV and MRP2. The study by Ray et al.43 showed that TFV was not substrate of MRP2 using transporter‐overexpressing cell lines. However, some gene association studies showed that ABCC2 variants were associated with TFV‐induced renal side effects in HIV‐infected patients.44, 45 One patient had very high TFV‐DP CLin in our studied cohort (408% of the median), potentially driving the significance of this relationship; however, the effects of our significant SNPs on this parameter were independent. FTC‐TP has higher CLout from the cells given an ABCC2_c.1249G>A variant, which increases MRP2 function.44, 45 The p16INK4a expression remains a significant predictor, potentially related to effects of cellular aging on the transport of drugs and/or the catabolic and anabolic processes for the metabolites.9

In the interaction analyses using linear regression models, we saw a significant effect of interaction between aging (both chronological age and p16INK4a expression) and gene polymorphism on ATV CLu/F and TFV‐DP CLin. However, we had two individuals with the ABCC4 variant in the ATV arm and two individuals with the SLC28A2 variant in the cohort, limiting the power of the current study. Other demographic variables were not included in the linear regression models, due to the small sample size. Hence, we performed Random Forest analysis that included all the available variables. Random Forest has been shown to have improved power over regression approaches to detect interactions in genetic association studies.46 However, due to the large number of predictors as compared with the limited sample size, the interaction effects of aging and polymorphisms may have been suppressed by other, stronger predictors in the Random Forest analysis.

Our study has some limitations. First, the small sample size, especially the ATV/RTV treatment arm, has limited our conclusions. We might have ruled out some important SNPs because of the low allele frequency. The low numbers of participants also compromised the distribution of some studied factors and the statistical power of our analysis; particularly, we lack participants with the frailty phenotype, limiting our ability to investigate frailty. Second, due to the large number of variables, we conducted the interaction analysis only involving aging (both chronological and cellular) and genetic factors, potentially excluding the interaction of other demographic factors with genetic factors. The interaction analyses were limited by low numbers of the minor alleles for the SNPs studied. The estimated effect sizes should be used to motivate well‐powered follow‐up studies to work to establish genetic factors as potential PK predictors.

Nevertheless, our study has several advantages. Using a novel statistical methodology, a population PK model is able to describe the PK processes more physiologically and explain the variability associated with those processes more appropriately than traditional PK analysis. Therefore, our analysis is based on more precise PK parameters, partially compensating for the small sample size. Additionally, we explored the SNPs at a gene‐wide level with the DMET array platform, maximizing the chance to identify relevant SNPs/genes. We found the associations of several SNPs with ARV PK that have not been reported before, generating interesting hypotheses for future study. Moreover, we evaluated the interaction between aging and SNPs in predicting ARV PK, using both the chronological and cellular aging. The significant interactions between aging and the genetic factors were seen for both chronological age and p16INK4a expression.

Our study provides data on the effect of aging on ARVs PK in the context of PGx, highlighting the importance of the cellular aging biomarker, p16INK4a, in predicting the ARV drug disposition and the potential alterations in aging‐PK relations by PGx. The investigations add knowledge to the current PGx of ARVs, and facilitate further studies with focus on precision treatments, especially those containing drugs with similar disposition profiles with our studied ARVs, for the growing older HIV‐infected population.

Conflict of Interest

The authors declared no conflict of interest.

Supporting Information

Click here for additional data file.

Supporting Information

Click here for additional data file.