Model-based tests to detect gene effect in pharmacokinetic studies
Bertrand J (1), Comets E (1), Laffont CM (2), Chenel M (3), Mentré F (1)
(1) UMR738, INSERM, Paris, France; Université Paris Diderot, Paris, France; (2) UMR181, Physiologie et Toxicologie Expérimentales INRA, ENVT, Toulouse, France; (3) Institut de Recherches Internationales Servier, Courbevoie, France
Introduction
Pharmacogenetics studies the relationship between the interindividual variability and variations in DNA sequence of proteins involved in mechanisms of drug absorption, distribution, elimination and effect [1]. Indeed, single nucleotide polymorphisms (SNPs) can change the amino acid sequence of proteins or else alter gene splicing, transcription factor binding, or the sequence of non-coding RNA. Pharmacogenetic information is more and more considered during the development process and the clinical use of a drug to realize the concept of personalized medicine. In this work we consider a biallelic SNP which leads to three possible genotypes all the more unbalanced when the mutant allele is rare.
In pharmacogenetic studies, concentrations data are mainly analysed using non-compartmental methods followed by a one-way analysis of variance (ANOVA) on the individual parameters of interest. More sophisticated approaches based on NonLinear Mixed Effects Models (NLMEM) have also sometimes been used with various approaches to test pharmacogenetic information. Preliminary screening is usually performed using ANOVA on the individual parameters estimates [2] followed by a stepwise model building approach with the likelihood ratio test (LRT) [3]. As an alternative approach, a global Wald test can assess whether estimates for the genetic effect are significant [4].
As in NLMEM the integral in the likelihood has no analytical form, specific algorithms are needed to estimate the model parameters and their standard error (SE). The two most widely used methods, first order (FO) and first order conditional estimation with interaction (FOCE-I), linearize the model function and are implemented in the NONMEM software [5]. The more recent stochastic EM algorithm (SAEM) avoids the linearisation step and is implemented in the MONOLIX software [6].
Aims
The aim of our work was to evaluate by simulation the three tests described above, ANOVA, LRT and Wald, in terms of type I error and power. We studied different estimation algorithms, we proposed type I error correction by means of permutation and we further investigated the impact of the design on the performances of these tests. We applied our conclusions to investigate the influence of SNPs on indinavir pharmacokinetics in HIV patients from the COPHAR2-ANRS 111 trial [7] and on concentrations of a drug under development (drug X) and its active metabolite.
Simulation study
a) Simulation settings
The concentrations were simulated using a one compartment model at steady state with first order absorption (ka), first order elimination (k), a diagonal matrix for the random effects and a proportional error model [8]. The model and parameters used for simulation came from a preliminary analysis without covariates of the indinavir concentration data from the COPHAR2-ANRS 111 trial. The genetic framework was inspired from the ABCB1 gene coding for the P glycoprotein found on the main physiological barriers. We simulated a diplotype of SNP1 (C and T the wild and mutant allele, respectively) and SNP2 affecting the drug bioavailability through the volume of distribution (V/F) The diplotype distribution mimicked that of exon 26 and exon 21 of the ABCB1 gene yielding for SNP1 unbalanced frequencies of 24%, 48% and 28% respectively for CC, CT and TT genotypes.
b) Influence of the estimation algorithm
Under the null hypothesis of no gene effect (H0), 1000 data sets were simulated with a design of N=40 patients and n=4 samples (inspired from the COPHAR2-ANRS 111 trial) and 1000 with the same sampling pattern but N=200 patients. Under the alternative hypothesis with a gene effect (H1), we simulated 1000 data sets with the N=40/n=4 design. Each of the three tests was applied to detect the effect of the SNP1 polymorphism. The type I error and the power of the tests were evaluated using FO, and FOCE-I implemented in NONMEM version V and SAEM implemented in MONOLIX version 2.1.
FOCE encountered many more convergence problems than FO whereas SAEM achieved convergence on all data sets. The ANOVA maintained a type I error close to the nominal level for both designs with all algorithms. The Wald test and the LRT obtained slightly inflated type I error on the N=40/n=4 design with both SAEM and FOCE-I. This inflation was corrected on the design N=200/n=4 for both tests with SAEM while only for the LRT with FOCE-I. The linearisation in the FO algorithm led to very poor performances of both the Wald test and the LRT.
c) Type I error correction
To correct type I error inflation of the Wald test and the LRT when N is small (40) we compared the approach by means of simulations and its non parametric alternative, the permutation test. We analyzed 200 data sets under H0 and 200 data sets under H1 for the design N=40/n=4. On the 200 data sets, the type I error and the power are estimated first using the theoretical threshold and then applying a threshold obtained from simulations or permutations. We used SAEM and FO but not FOCE-I to limit computing times and numerical difficulties.
Both correction approaches provided a type I error close to the nominal level with FO or SAEM. However, the corrected power for the Wald test and the LRT were much lower using FO than with SAEM. Indeed, we observed with FO a high correlation between the estimation errors of the parameter of the gene effect and their estimates, this relationship led to decreased values of the Wald statistic and therefore reduced the power to detect a genetic polymorphism effect. This pattern was not observed using SAEM. Permutations of the genotypes vector could not overcome limitations of the estimation method.
d) Influence of the design
In 2007, the EMEA has stated that pharmacogenetic studies should include a satisfactory number of patients of each geno- or phenotype in order to obtain valid correlation data [9]. Therefore, in a third step, we considered two other designs with a larger number of subjects but different blood sampling strategies as four sampling on each patient would no longer be practical. A design optimized using the PFIM software [10], N=80/n=2 sorted in four groups and a combined design, N=20/n=4 plus N=80 with only a trough concentration. These two designs involved the same total number of observations as the original design with 40 patients, to allow proper comparisons between designs [11]. In front of the results obtained previously, the estimation was performed only with the SAEM algorithm.
ANOVA kept a correct type I error estimate but its power was lowered when there was high shrinkage (N=100/n=4,1). The Wald test and the LRT had significantly increased type I error in the two designs. Yet, the inflation remained moderate as all the type I errors were below 10%. The design N=80/n=2 provided the best performances; as it had both the lowest estimation errors on the gene effect coefficients and the highest power among the three designs with a total of 160 observations.
Applications
a) Influence of five SNPs on indinavir pharmacokinetics in the COPHAR2-ANRS 111 trial
We analysed the concentrations of forty patients collected after two weeks of treatment at 1, 3, 6 and 12 hours following administration. For each patient, genotypes are obtained for ABCB1 exon 26 and 21, CYP3A5*3, CYP3A5*6 and CYP3A4*1B polymorphisms along with demographic covariates. We modelled the pharmacokinetics of indinavir using SAEM in MONOLIX version 2.1. We performed a screening on the EBE using ANOVA followed by an ascending selection based on LRT. In front of the previous results, as the number of patients is limited, we then performed a backward selection using a LRT by permutation [12].
A one-compartment model with first-order absorption and elimination best described the indinavir concentrations. The eight patients *1B/*1B for CYP3A4 gene had an absorption rate decreased by 70% compared to *1A/*1B or *1A/*1A genotypes (0.5 versus 2.1, P=0.04).
b) Influence of the CYP2D6 polymorphism on the pharmacokinetics of the drug X and its active metabolite
We analysed the pharmacokinetic profiles of ninety patients collected after four weeks of treatment at 1, 3, 6 and 24 hours following administration. Concentrations were measured for both the parent drug and its active metabolite. For each patient, genotypes are obtained for the CYP2D6 polymorphism. We modelled jointly the concentrations of drug X and its active metabolite using SAEM in MONOLIX version 2.4. We performed a screening on the EBE followed by an ascending selection based on LRT. As the design of the study was rather rich, we did not perform a correction using permutations.
In the final model, the volumes of the parent drug and the metabolite were set to be equal and the decrease of bioavailability with the dose was taken into account. Through a first pass effect, 14% of the dose formed directly the metabolite and this metabolite was partly back-transformed into the parent drug. Between patients variabilities were estimated on all parameters with the exception of the metabolization and back-transformation clearances while correlations were observed between the absorption constant rate, the volume of distribution and the elimination clearance of the metabolite. The CYP2D6 polymorphism was found to influence the pharmacokinetics of the molecules.
Conclusions
For test of gene effects in pharmacokinetic studies analysed with NLMEM, we recommend using an ANOVA on the individual parameters as this method showed, in our simulation setting, the best performance in terms of type I error whatever the estimation method even in designs with some shrinkage. If one still wants to perform Wald test or LRT and construct a global model, we promote the use of permutation tests on designs with unbalanced genotypes and/or small number of subjects. A more recent estimation method, as the SAEM algorithm, shows better properties and is easier to apply for permutation tests because it ensures reasonable computing time and meets with no numerical difficulties on repeated data sets.
References
[1] EMEA. ICH topic E15 definitions for genomic biomarkers, pharmacogenomics, pharmacogenetics, genomic data and sample coding categories. Technical report, EMEA, 2008.
[2] D Hirt, F Mentré, A Tran, E Rey, S Auleley, D Salmon, X Duval, JM Tréluyer and the COPHAR2-ANRS Study Group. Effect of CYP2C19 polymorphism on nelfinavir to M8 biotransformation in HIV patients. Br. J. Clin. Pharmacol., 65:548-57, 2008.
[3] Y Yamasaki, I Ieiri, H Kusuhara, T Sasaki, M Kimura, H Tabuchi, Y Ando, S Irie, JA Ware, Y Nakai, S Higuchi, and Y Sugiyama. Pharmacogenetic characterization of sulfasalazine disposition based on NAT2 and ABCG2 (Bcrp) gene polymorphisms in humans. Clin. Pharmacol. Ther., 84:95-103, 2008.
[4] D Li, W Lu, JY Zhu, J Gao, YQ Lou, and GL Zhang. Population pharmacokinetics of tacrolimus and CYP3A5, MDR1 and IL-10 polymorphisms in adult liver transplant patients. Clin. Pharm. Ther., 32:505-515, 2007.
[5] L Sheiner and S Beal. NONMEM Version 5.1. University of California, NONMEM Project Group, San Francisco, 1998.
[6] M Lavielle. MONOLIX (MOdèles NOn LInéaires à effets miXtes). MONOLIX group, Orsay, France, 2008. http://software.monolix.org/index.php.
[7] X Duval, F Mentré, E Rey, S Auleley, G Peytavin, M Biour, A Métro, C Goujard, AM Taburet, C Lascoux, X Panhard, JM Tréluyer, D Salmon. Benefit of therapeutic drug monitoring of protease inhibitors in HIV-infected patients Depends on PI used in HAART Regimen - ANRS 111 trial. Fund. Clin. Pharmacol., in press, 2009.
[8] J Bertrand, E Comets, and F Mentré. Comparison of model-based tests and selection strategies to detect genetic polymorphisms influencing pharmacokinetic parameters. J. Biopharm. Stat., 18:1084-1102, 2008.
[9] EMEA. Reflection paper on the use of pharmacogenetics in the pharmacokinetic evaluation of medicinal products. Technical report, EMEA, 2007.
[10] S Retout, E Comets, A Samson, and F Mentré. Design in nonlinear mixed effects models: optimization using the Fedorov-Wynn algorithm and power of the Wald test for binary covariates. Stat. Med., 26:5162-5179, 2007. http://www.pfim.biostat.fr/.
[11] J Bertrand, E Comets, and F Mentré. Properties of different tests to detect the effect of a genetic covariate on pharmacokinetic parameters using the SAEM algorithm for several designs. PAGE 17 [Abstr 1337], Marseille, France, 2008.
[12] J Bertrand, JM Treluyer, X Panhard, A Tran, S, E Rey, D Salmon-Céron, X Duval, F Mentré and the COPHAR2-ANRS 111 study group. Influence of pharmacogenetics on pharmacokinetic interindividual variability of indinavir and lopinavir in HIV patients (COPHAR2 - ANRS 111 trial). PAGE 16 [Abstr 1153], Copenhagen, Denmark, 2007.
Acknowledgments
During this work, Céline M. Laffont was working at the Institut de Recherches Internationales Servier (France) as pharmacometrician and Julie Bertrand was supported by a grant from the Institut de Recherches Internationales Servier (France).
We would like to thank the COPHAR 2-ANRS 111 scientific committee (investigators: Pr. D. Salmon and Dr X. Duval, pharmacology: Pr JM. Tréluyer, methodology: Pr F. Mentré) for giving us access to the pharmacogenetic data of the indinavir arm. We would also like to thank the IFR02 of INSERM and Hervé Le Nagard for the use of the "centre de biomodélisation'' as well as the Pr. Marc Lavielle for the precious help he provided in using MONOLIX.
