Pharmacokinetic Based Design of Individualized ... - ACS Publications

ARTICLE pubs.acs.org/IECR

Pharmacokinetic Based Design of Individualized Dosage Regimens Using a Bayesian Approach Jose Miguel Laínez,*,† Gary Blau,† Linas Mockus,† Seza Orc-un,‡ and Gintaras V. Reklaitis† † ‡

School of Chemical Engineering, Purdue University, West Lafayette, Indiana 47907, United States e-Enterprise Center, Discovery Park, Purdue University, West Lafayette, Indiana 47907, United States ABSTRACT: Clinical trials and health care studies generate an enormous amount of data. This data is used by pharmaceutical companies during new drug development processes to characterize patient populations and determine a standardized dosage regimen for new patients, make commercial decisions, and gain approval from regulatory agencies. Nevertheless, the knowledge embedded in such data is rarely further exploited for customized patient care. In most cases, there is a significant difference between the pharmacokinetic profile of patients in a population, yet these differences are not reflected in the standardized dosage regimen. Here, a Bayesian methodology is proposed to individualize dosage regimens by combining the pharmacokinetic data collected from a patient population during clinical trials and additional data coming from a minimal number of serum samples from the new patient. In the Bayesian sense, the distribution of pharmacokinetic parameters from the population data is treated as prior information, and the posterior patient specific distribution of pharmacokinetic parameters is calculated. Then, such a posterior distribution is used to obtain dosage regimens that result in drug concentrations that are kept within the therapeutic window at a target confidence level for that patient. Moreover, a methodology is presented to suggest the sampling schedule for new patients so as to reduce the number of samples required to obtain well characterized individual pharmacometric parameters. Available pharmacokinetic data for Gabapentin, a therapeutic agent for epilepsy and neuropathic pain, is used to illustrate the concepts underlying the proposed strategy and the benefits of an individualized regimen over a standardized dosage regimen.

’ INTRODUCTION Every pharmaceutical company must submit a new drug application to the Food and Drug Administration (FDA) to obtain marketing approval for a drug within the U.S. Similar drug approval processes are implemented in other regions and countries as well. The new drug application contains data regarding the safety and efficacy of the drug which was obtained during the research and development phase. This data includes the results of clinical trials, pharmacology and toxicology data, chemistry and manufacturing data, and proposed packaging and labeling information.1 The pharmacometric information derived to obtain marketing approval is mostly used to determine a “one dose fits all” or standardized dosage regimen (i.e., dose amount and interval of administration) which is often sought by the drug sponsor because it results in simplified production/packaging processes and ease of use by practitioners.2 However, given that individuals vary significantly in their response to drugs, a “one dose fits all” approach may result in ineffective or excessive individual dosing and additional costs due to overmedication. The latter point is especially important for those drugs with narrow therapeutic windows and highly toxic effects. In fact, a substantial proportion of patients do not respond at all, respond only partially, or experience adverse drug reactions for all major classes of drugs given at standard doses.3 Problems associated with drug therapy are frequent and costly. According to a recent survey by Krahenbuhl et al.,4 adverse drug events affect 6.1% of hospitalized patients. Among the wards in which the patients were studied, the frequency was higher in geriatric wards, with a median of 27.9%. Indeed, one of the main factors associated with therapy problems is age: Steed5 suggests that in Canada overmedicated seniors are clogging emergency departments for no apparent reason. r 2011 American Chemical Society

Starting with the standardized dosage regimen supplied by the manufacturer, a physician typically has to resort to titration studies or dose ranging trials to determine a suitable dosage regimen for an individual patient. A titration study, a protocol consisting of a combination of experiments and decision rules, constrained by toxicities is conducted as follows: if the previous observation indicates toxicities, then the dose is decreased by one predetermined increment. If no toxicities are observed, the dose is increased by one predetermined increment.6 Usually, the patient needs to visit the physician for dose adjustments until the desired therapeutic effect is achieved. Not only is this process expensive and cumbersome, but it is also easy to overmedicate a patient, particularly if the side effects are less pronounced. A different framework for selecting the components of an individualized dosage regimen is provided by a pharmacometric approach, in which decision making relies on the application of mathematical and statistical methods. Pharmacometrics can be broadly divided into two major areas, pharmacokinetics (PK) and pharmacodynamics (PD). The former studies how the body handles a drug by describing the time course of drug concentration in plasma, while the latter describes the relationship between drug concentration and therapeutic effect. The utility of this approach is that it recognizes the importance of relating drug effects to plasma concentrations, rather than exclusively to dosage.7 Special Issue: Puigjaner Issue Received: July 28, 2010 Accepted: January 24, 2011 Revised: January 3, 2011 Published: March 04, 2011 5114

dx.doi.org/10.1021/ie101610r | Ind. Eng. Chem. Res. 2011, 50, 5114–5130

Industrial & Engineering Chemistry Research PK models assume that the major mechanisms influencing the dynamic deposition of a drug in the body are its absorption, distribution, metabolism, transfer across plasma membranes, and elimination.8 In its simplest linear form, where the various anatomical components of the body are assumed to be well mixed compartments, the mechanism may be characterized by a mathematical model represented as a sum of exponential terms as stated in eq 1. In this equation, the set of PK parameters is comprised of parameters As and Bs, and Ct represents the predicted drug concentration at time t. X ^t ¼ C As expð - Bs tÞ ð1Þ s

The linearity refers to the fact that the drug concentration change with respect to time is directly proportional to the concentration. Nonlinear models usually require the solution of differential algebraic equations to predict the concentration at a specific point in time starting from some initial conditions. An example is the Michaelis-Menten rate model.9 On the other hand, the basic idea of physiological PK is to extend compartmental PK modeling so that quantitative aspects of other biological areas can be incorporated. This approach focuses on models that anatomically describe real local tissue regions, including their blood flow, binding, and transport characteristics. The usual criticism of the use of compartmental models is that the estimated parameters are meaningless. However, Wagner9 has demonstrated that there exists a relationship between physiological and compartment based models. While this is possible, it is not generally the case. The relationship between physiological and compartment based models is often discovered with posthoc correlations. In this paper, we will not cover the methodology for selecting the most suitable PK model that best describes the concentration time data for a treatment agent. Blau et al.10 have addressed this problem using Bayesian methods. Instead, we will assume that the model is known and that the parameters of the PK model are different between patients so that individualizing the dosage regimen makes a difference. Therefore, we will deal with only the following issues: (i) determining the prior distribution of the parameters for the population from the data collected during drug development where the model is known; (ii) determining the location and number of data points needed to provide an adequate estimate of the individualized patient plasma concentration time profile; (iii) calculating the posterior distribution for an individual using this prior distribution and plasma concentration time data collected on that individual; (iv) predicting the concentration time profile for an individual using the standardized dosage; and (v) determining a dosage regimen which maximizes the efficacy while minimizing the toxicity when the therapeutic window is available. The above process will be illustrated in a retrospective study on an important treatment agent, Gabapentin. The paper is organized as follows. The section entitled Population Pharmacokinetics formally describes the problem of population PK and presents a literature review of the utilization of population PK and dosage regimen individualization. In The Proposed Bayesian Framework section, the general framework of the Bayesian based methodology is explained. Here, the building of the Population Prior is presented as well as the procedure to determine the sampling schedule required to get well characterized individual PK parameters. In the section Stage III: Dosage Regimen Optimization, it is shown how the patient parameter distribution can be exploited to optimize the dosage regimen design. The

ARTICLE

performance of the proposed methodology is illustrated in the section Case Study: Gabapentin. Finally, some conclusions about this work are drawn in the Final Remarks section.

’ POPULATION PHARMACOKINETICS When formulating a population model, it is recognized that it is common to find large interindividual variation in the PK of a drug so that the overall variability in the experimental measurements/observations reflects both measurement error and interindividual variability.11 Then, a PK model within the framework of a population can be described as follows: 9 ^yj ¼ fj ðxj ðtÞÞ > > > > > dxj ðtÞ = ¼ gj ðxj ðtÞ, uj , φÞ "j ∈ f1; :::; J g ð2Þ dt > hj ðxj ðtÞ, uj , φÞ ¼ 0 > > > > ; x ðt Þ ¼ x j

o

oj

where ^yj is an I-dimensional vector of predicted values for the observations for the jth patient (e.g., drug blood concentration, urinary drug levels). In this model fj is a vector function relating the state variables xj(t) to the predicted observations for patient j. The vector φ represents the unknown PK parameters, and uj is the vector of controlled variables for patient j. The controlled variables are set in the clinical trials, and their numerical values are assumed to be precisely known. The vectors of differential and algebraic equations that define the PK model are gj and hj, respectively. Finally, xoj is the vector of initial conditions. Different doses, administration schedules (e.g., single, multiple doses), dosage regimens, and routes of administration are commonly used during clinical trials. In addition, several kinds of responses (e.g., blood and urinary drug levels) might be measured. Accordingly, the model could differ among patients. Nevertheless, it is assumed that the set of underlying parameters is qualitatively the same among patients.11,12 yj ¼ ^yj þ εðxj , uj , ϑÞ "j ð3Þ Assuming additive error, the experimental observations for the jth patient can be expressed as in eq 3. Here, ε is the random vector representing the experimental error in the measured values, and ϑ is the vector of parameters defining the experimental error variability. The Experimental Error. It will be assumed that measurement errors are independent and normally distributed. ð4Þ N ðεj0, σ 2 ðxj , uj , ϑÞÞ In order to accommodate heteroscedasticity in the measurement errors inherent to concentration measurements, Reilly and Blau13 have proposed the following error model: σ2 ðxj , uj , ϑÞ ¼ ω2 ð^yj Þγ þ ζ2

ð5Þ

where ω, γ, and ζ are the parameters defining the error model. Homoscedastic errors can be accommodated by this model by letting γ take a value equal to zero. In the particular case of PK models where the drug is administered orally and rapidly absorbed into the bloodstream, measurement errors are exacerbated if the sample time is not recorded precisely. To deal with this, error variance models which incorporate time dependency have been proposed in the 5115

dx.doi.org/10.1021/ie101610r |Ind. Eng. Chem. Res. 2011, 50, 5114–5130

Industrial & Engineering Chemistry Research

ARTICLE

works of Ruppert et al.14 and Chan et al.15 This will not be considered in our development. Population Pharmacokinetic Parameters. Assuming that the PK model represented by eq 2 is valid and the set of underlying parameters is the same among patients, the PK parameters for any individual j can be written as φj ¼ φh þ ηj

population mean. However, the individuality of each subject is maintained. The most widely used approach is the first order method which linearizes the model in the random effects (i.e., experimental error and intersubject variability) by using a first order Taylor series expansion. The parameters are then commonly estimated using a maximum likelihood approach.20 Unfortunately, linearization introduces a severe limitation on the method, resulting in biased parameter estimates with unrealistic variances. Bayesian methods overcome this limitation and permit the PK and error models to exhibit their nonlinear behavior. Nonparametric Maximum Likelihood. Mallet21 proposed a method which did not assume that the population parameters distribution belonged to a given parametric family of distributions. A maximum likelihood approach is used to estimate the multivariate discrete distribution of the parameters of a regression model among a population, from individuals’ samples. The method assumes the experimental error distribution is known. Mallet confirmed the ability of this method to reveal unexpected features in the population parameters distribution such as multimodalities by using a simulated example. Seminonparametric Method Approach. This method is based on an assumption that the distribution of the random effects has a smooth density (i.e., it has continuous derivatives up to some desired order over some domain). The smooth distribution can be represented as an infinite series expansion. This approach uses a finite number of terms resulting from an approximation of the infinite expansion. The number of terms to be used is determined by a single tuning parameter.22

ð6Þ

where φh denotes the vector of mean values of the PK parameters for the population and ηj is the variation of the individual j parameters around the population mean values. Then, the individual parameters are assumed to arise from a multivariate probability function P as shown in eq 7: P ðφÞ ¼ F ðφjφh, ΩÞ

ð7Þ

Here, F can be characterized by φh, the vector of mean population PK parameters, and Ω, the variance-covariance matrix of between-subject random variability. It is also often useful to consider the values of the parameters for the individual to be related to the population parameter by means of covariates, in which case the expression may be defined as expressed in eq 8. In this equation, w denotes the relationship between the covariates zj (e.g., age, weight) and the vector of the mean population PK parameters.16 P ðφÞ ¼ F ðφjwðφh, zj Þ, ΩÞ

ð8Þ

Different methods have been proposed to determine the population PK parameters. These methods can be categorized as follows:12,17,18 The Naive Average Data Approach. This method requires identical sampling schedules for all patients. The average value of observations at each sampling time is computed. Then, the population PK parameters are determined by fitting the PK model to the mean-observations vector. The Naive Pooled Data Analysis. This method assumes that all observations arise from one unique patient. This reference patient is characterized by a set of PK parameters that is considered the population PK parameters. Two Stage Approach (TSA). The first stage of this method consists in determining the PK parameters for each individual by fitting the PK model to the observations corresponding to each individual, provided there are sufficient observations from an individual. In the second stage, φh and Ω are determined on the basis of the individual PK parameters. There are various methods which differ in the procedure implemented in the second stage. For instance: • The standard TSA estimates φh as the mean of the individual’s PK parameters. In this approach, Ω is considered to be a diagonal matrix whose values are the linearized variance of a patient’s PK parameters around the previously calculated φh. • A more sophisticated method is the global TSA. This approach assumes that the individual parameters are normally distributed around the mean of the population PK parameters. Then, φh and Ω are determined following a maximum likelihood approach.19 Nonlinear Mixed-Effects Modeling (NONMEM). This methodology analyzes all data from all subjects together with patient specific information. This method uses a hierarchical model to directly estimate the parameters of the population from the full set of observations. The individual’s parameters are assumed to be normally distributed and vary around the

The drawback of the naive average data approach and the naive pooled data analysis is that they ignore the interpatient variability. They do not distinguish between individual differences and other sources of variability. The drawback of TSA and the nonlinear mixed-effects modeling approaches is that they make assumptions with regard to the shape of the population parameter distribution. This tends to overestimate PK parameters dispersion, which is associated with the biases in the variance-covariance matrix Ω.11 Although the nonparametric likelihood approach does not require the assumption of a shape for the distribution of population parameters, it uses maximum likelihood point estimates for each patient to build a discrete parameter distribution for the population. Finally, the shape of the distribution obtained by using the seminonparametric approach is restricted by the tuning parameter, which is selected a priori. In the section The Population Prior, we describe a methodology to overcome these disadvantages in the context of a Bayesian approach. Previous Work. The potential impact that PK modeling can have on finding optimal individualized treatments has been widely recognized. In the mid-1980s, Whiting et al.23 pointed out that the dosage of drugs chosen on a relatively empirical basis may not lead to optimum treatment. They proposed to adjust the dose given to patients by analyzing early plasma samples. A NONMEM approach was used to estimate the patient PK. They indicated that these estimates allowed selection of future doses which were more appropriate for each patient. Unfortunately, this procedure was not described. Jelliffe et al.18 emphasized the importance of population PK models in order to store past experience of drug behavior. Then, they devised estimates of an individual’s parameters by minimizing the combination of deviations from population parameters and predicted plasma concentrations to make 5116



ARTICLE

Figure 1. Bayesian framework for individualized dosage regimen optimization.

an individualized patient-specific model. They suggested that adaptive control could be used to develop a dosage regimen which would achieve the desired therapeutic goal. Mehvar24 showed how PK principles could be used to design dosage regimens for individual patients in order to achieve therapeutic plasma concentrations of drugs. He did not follow a population approach; instead the PK parameters were fitted to plasma samples taken from the patient. Then, a dosage regimen was determined so that the average of plasma concentration at a steady state was equal to an estimated target. In the spirit of the Bayesian approach that we advocate in this paper, Duffull et al.25 compared the population model obtained by using a hierarchical Bayesian, described in the next section, and a NONMEM approach. They demonstrated the superiority of the former over the latter in capturing the between-subject variability. However, they did not address the dose individualization problem. A procedure to determine individualized dosage regimens was presented by Salinger and co-workers.26,27 This procedure intended to define a dosage regimen which achieved a target area under the concentration-time curve by using a maximum posterior distribution fitting.28 Barrett et al.29 developed a visualization tool to support drug-specific decisions related to the pharmatherapeutic management of individual patients. The system used a NONMEM approach with its limitations to generate individual patients’ profiles. Salinger et al.30 used a maximum posterior distribution fitting to determine an outpatient sampling schedule for intravenous Busulfan administration. The optimal schedule is selected from a potential pool of sampling times on the basis of a scaled mean squared error of the PK parameters. Most of the procedures to date for determining an individualized dosage regimen are based on the use of point estimates rather than probability distributions. Usually, a target is defined (e.g., expected area under the concentration-time curve, steady state average concentration), which must be achieved by the

optimal individual dosage regimen. By contrast, the work described in the following sections proposes a dosage regimen optimization procedure which relies on probability distributions of the parameters and predicted concentrations. Such distributions are the result of following a fully Bayesian methodology. The dosage regimen generates concentration distributions which are within the therapeutic window with a target confidence level. Additionally, we propose a methodology for determining a minimum number of samples and their timing (relative to the drug administration) required in order to characterize the joint distribution of PK parameters for a new patient.

’ THE PROPOSED BAYESIAN FRAMEWORK The Bayesian approach presented in this work follows the general mathematical model building framework described by Blau et al.10 When a frequentist interpretation of statistics is used in model building, it is considered that there is only one single true set of parameters that fits the experimental data. Only that data is used to select a model and determine parameter estimates (e.g., regression analysis, maximum likelihood approaches). In other words, the data “speaks for itself” and is unbiased by the experimentalist or model builder. In the Bayesian approach, however, prior knowledge available to these individuals is used to interpret results. The use of prior knowledge is of interest to us since we intend to take advantage of the “prior” population data that is generated during clinical trials to obtain the distribution of PK parameters specific to an individual from only a few observed data points. Using a Bayesian analysis, it is possible to generate posterior joint probability distributions for the PK parameters. In this work, we assume that an adequate PK model is selected during the drug development process and focus on generating tighter probability distributions of the kinetic parameters for new 5117



ARTICLE

patients to drive individualized dosage regimens. As alluded to previously, model selection falls nicely into a Bayesian framework but is beyond the scope of this article. We refer the interested reader to Blau et al.10 and Hsu et al.31 The framework we propose for individualized dosage regimen optimization entails three stages, as depicted in Figure 1. The first stage comprises the use of clinical trials data for determining the joint posterior distribution of PK parameters for each individual involved in such clinical studies using a Bayesian approach. At this point, the inference is just driven by the data collected during clinical trials. It is assumed that the clinical trials data is of sufficient quality to yield a good characterization of these individuals. Then, these individual distributions are used to construct a joint probability distribution of the PK model parameters for all of the individuals, which will be called the population prior. A methodology has been developed to produce a sampling schedule (i.e., number of samples and their corresponding collection times) to obtain well characterized PK parameters for a new patient. This methodology uses as input the population prior. This first stage can be regarded as an “offline” process since it is performed just once to formalize the prior knowledge required by the Bayesian approach. The second and third stages are repeated for each new patient for whom the dosage regimen is to be individualized. In the second stage, the posterior joint distribution of PK parameters is computed for the new subject using the population prior and the data collected for the patient at the times following administration suggested in stage I. In the third stage, the posterior joint distribution of PK parameters is used to compute a dosage regimen that places the plasma concentration within the therapeutic window given a target confidence level. For this stage, the dosage regimen design is posed as a constrained optimization problem. It will be shown that by fixing the time between doses the problem is dramatically simplified. This is a suitable approach for determining a dosage regimen whose dosing interval can be easily followed by the patient. Next, we describe each stage in detail. Stage I. Assume that, to initiate this process, data collected on a population or cohort of patients (Dj) and a PK model validated with this data are needed. Parameter Estimation. The PK parameters are computed for each patient j involved in the clinical trials. In accordance with Bayes’ theorem (eq 9), the posterior beliefs about the PK parameters of patient j involve the incorporation of two elements, namely, (i) prior beliefs, given by the prior distribution p(φ,ϑ), and (ii) the experimental outcomes. The probability of achieving the experimental outcomes given the model and a given set of parameters is captured by the likelihood function (L(Dj|φ,ϑ)). LðDj jφ, ϑÞ pðφ, ϑÞ ð9Þ φ ϑ LðDj jφ, ϑÞ pðφ, ϑÞ dφ dϑ As discussed previously, it is assumed that the experimental errors for each of the i available data points for patient j are independent and normally distributed with mean 0 and variance σi2. Then, the probability for the data points included in the data set Dj can be expressed as in eq 10. Y pðεij Þ pðεj Þ ¼ pj ðφ, ϑjDj Þ ¼ R R

i ∈ f1:::Ig

¼

Y i ∈ f1:::Ig

( 1 ð2πÞ1=2 σ i

exp -

ε2ij 2σ 2i

!) "j

ð10Þ

According to eq 3, for a given set of PK parameters, the experimental errors are given by εij ¼ yij - ^yij "i, j

ð11Þ

Then, the likelihood function can be expressed as follows: 8 Y

= L TW TW U e Dose e max Dosemin , L > min Q φ, R ðt, τÞ> max Q Uφ, R ðt, τÞ : ; t ∈ ½0, τ t ∈ ½0, τ

ð25Þ

Unfortunately, similar simplifications cannot be obtained for the case in which the interval is optimized for a given dose amount. A Global Approach. Instead of attaining the confidence level R at every time t ∈ [0,τ] as stated in eq 18, we may be interested

φ ðτÞ < Q U ðτÞÞ ¼ R L ðτÞ < Q PrðQ φ, R φ, R

ð26Þ

φ is given by the Here, the density function/probability of Q probability of obtaining a given value of Q φ at any time point within the dosing interval τ as shown in eq 27. Z τ pðQ φ ðt, τÞÞ dt φ ðτÞÞ ¼ 0 ð27Þ pðQ τ The dosage regimen optimization problem for the global approach can be then posed as follows: U ðτÞ - Q L ðτÞÞ min DoseðQ φ, R φ, R

Dose, τ

L ðτÞ g TW L Dose Q φ, R Dose QφU, R ðτÞ e TW U Dose g Dmin τ g τmin

ð28Þ

In the case of a specified interval τ, a feasible interval for the dose amount is given by eq 29. Notice that for this global approach, L U φ,R φ,R Q (τ) and Q (τ) each represent a single point rather than a function. ( ) L TW TW U ð29Þ max Dosemin , L e Dose e U ðτÞ ðτÞ Q Q φ, R φ, R

’ CASE STUDY: GABAPENTIN Clinical trial data from a study of Gabapentin is used to illustrate the proposed methodology. This data has been collected during the clinical studies carried out by Urban et al.35 Gabapentin is an anticonvulsant that is prescribed for epilepsy and other neuropathic disorders. The data available is for 19 patients who participated in a 36 h study of Gabapentin PK. A single dose of Gabapentin (400 mg) was administered orally. Blood samples were collected at 14 time points after dosing (0.25, 0.75, 1.00, 1.50, 2.00, 3.00, 4.00, 6.00, 8.00, 12.00, 18.00, 24.00, 30.00, and 36.00 h). The Model. The PK model selected for predicting plasma concentration is the single dose one-compartment model with first order absorption. The equation for this model is as follows: ka ^ t ¼ Dose F ½expð-ke ðt - to ÞÞ - expð-ka ðt - to ÞÞ C V ka - ke ð30Þ Here, F is the bioavailability (i.e., the fraction of the dose that reaches the circulation), V is the volume of distribution, Dose is the administered dose amount, ke is the elimination rate, ka is the absorption rate constant, and to is the dead time. In eq 30, t is the time after dosing or administration. For this case study, the input variables are t and Dose. Notice that the predicted concentration ^ t should be equal to zero for all t < to. The set of PK parameters C (φ) is given by (F/V), ke, ka, and to. F/V is used as a PK parameter since F is unknown. By doing so, we avoid having a so-called flipflop model (i.e., a model whose predictions are identical for two or more sets of solutions of the parameters). Notice that if clearance is used as a parameter, F and V must be considered as 5121



ARTICLE

Table 1. Mean and Standard Deviation (SD) of PK Parameters for Each Patient F/V (10-3/l)

ke (h-1) mean

SD

ka (h-1) mean

SD

to (h)

patient

mean

SD

mean

SD

P01

11.704

2.286

0.232 0.064 0.788 0.210 0.265 0.059

P02

10.047

0.552

0.194 0.019 1.469 0.171 0.241 0.022

P03 P04

13.209 9.297

0.995 1.619

0.207 0.023 0.871 0.101 0.244 0.026 0.129 0.030 0.498 0.127 0.208 0.060

P05

7.577

1.286

0.140 0.035 0.723 0.191 0.198 0.055

P06

7.799

1.486

0.193 0.070 1.694 0.636 0.232 0.055

P07

9.815

0.867

0.175 0.024 0.860 0.115 0.261 0.033

P08

9.493

1.755

0.199 0.050 0.709 0.190 0.182 0.049

P09

5.359

0.134

0.133 0.008 4.277 1.635 0.245 0.008

P10

12.560

1.987

0.255 0.050 0.586 0.123 0.171 0.039

P11 P12

13.766 12.540

2.292 2.681

0.176 0.036 0.329 0.075 0.298 0.087 0.229 0.069 0.732 0.219 0.285 0.073

P13

12.955

0.850

0.171 0.017 0.730 0.080 0.134 0.028

P14

9.379

1.371

0.121 0.030 0.966 0.304 0.190 0.053

P15

11.416

1.923

0.183 0.039 0.422 0.092 0.294 0.063

P16

18.495

2.757

0.221 0.040 0.427 0.081 0.265 0.051

P17

7.752

0.537

0.115 0.013 0.685 0.086 0.164 0.039

P18

10.733

1.805

0.250 0.061 0.971 0.220 0.269 0.045

P19

9.258

1.240

0.169 0.032 0.748 0.152 0.156 0.040

separated parameters in order to model the elimination process. In cases where covariates and their relationship with clearance and V are accounted, additional data to estimate the bioavailability (F) must be considered to avoid a flip-flop model. With regard to the error, there are no replicates which could be used to define an appropriate model for it. Therefore, a constant variance has been assumed for all data points (i.e., homoscedasticity). Although, an adequate model is assumed as an input to our methodology, the model adequacy can be confirmed by using the global lack of fit test described by Blau et al.10 This test compares the occurrence of the experimental points within highest probability density (HPD) regions for predictions (i.e., distribution of blood concentrations in our case). Using a confidence level of 95% for the HPD concentration region, 262 out of 265 experiments are found to lie within it, which results in a confidence level for the lack of fit of 0.014%. Hence, the selected model is adequate as measured by the global lack of fit test. The Population Distribution of PK Parameters. The methodology starts by computing the PK parameters for each of the individuals for which clinical trial data exits, by applying the Bayesian approach. A log-normal noninformative prior distribution has been assumed for the PK parameters (see the Parameter Estimation section). As previously mentioned, a constant variance has been assumed for the drug concentration. Table 1 summarizes the mean and standard deviation for each PK parameter and patient. It is evident that there exists intersubject variability. For Gabapentin, individuals vary more significantly in their absorption rate (ka) and bioavailability over volume of distribution (F/V). The individual’s distributions of PK parameters are used to build the population prior P as explained in The Population Prior section. The population prior has been built by (i) using the PK parameters samples for each individual and (ii) approximating each individual PK parameter’s distribution using a

Figure 4. Comparison of population prior approximation. Marginal distribution of PK parameters (dashed line, approximation by using a mixture of normal distributions; solid line, empirical distribution).

multivariate normal distribution. Figure 4 shows the population prior marginal distributions for both cases. As can be seen, the mixture of multivariate normal distributions renders a good approximation which will be used for further computations. Sampling Schedule Selection. Here, the methodology previously described for determining a sampling schedule is followed. Once, the population prior P has been built, it is used to obtain the distribution of blood concentration at different time points. It is assumed that the PK study for new patients is carried out by administering a single dose. We have selected 0.75, 1.00, 1.25, 1.50, 1.75, 2.00, 2.50, 3.00, 3.50, 4.00, 6.00, 8.00, 10.00, 12.00, 14.00, 18.00, 24.00, and 30.00 h after drug administration as potential time points for blood sampling. Earlier sample points are preferred over later points since they would allow a shorter stay for the patient in the clinic. Figure 5 depicts the predicted blood concentration distributions for each considered time point by using the population prior. Table 2 summarizes some variability measures for the most relevant time points. Blood samples at 1.50, 2.00, and 2.50 h after dosing are good candidates for the first blood draw as the drug concentration at those times presents higher variability. Consequently, these three times are tested as candidates for the first sampling point. The predicted drug concentrations for these times are simulated by using the maximum density estimates of PK parameters from the population prior (P ) which are F/V = 7.667 10-6 mL-1, ke = 0.115 h-1, ka = 0.656 h-1, and t0 = 0.121 h. The simulated drug concentrations for 1.50 h, 2.00 h, and 2.50 h are 1668.6, 1911.8, and 2047.2 μg/mL, respectively. After applying the Bayesian approach to obtain a posterior distribution for each of the simulated data by taking P as prior, the time point that reduces the predicted drug concentrations’ variability the most is selected. The variability of the concentrations obtained by considering each of these time points as the first sample is summarized in Table 3. Note that the sample at time 1.50 h turns out to be the most appropriate point to select (it results in the lower average standard deviation, 195.5 ng/mL). The probability distributions of drug concentration obtained by using the sample at time 1.50 h and the population prior (P ) are shown in Figure 6. Notice that the concentration variability has been significantly reduced by only considering one data point, the drug concentration at 1.50 h after dosing. However, important improvements 5122



ARTICLE

Figure 5. Predicted drug concentration (ng/mL) distribution for different time points by using the population prior.

Table 2. Variability Summary for the Most Relevant Points by Using the Population Prior time after dosing

standard deviation (ng/mL)

mean (ng/mL)

CV

t = 1.50

464.537

2170.372

0.214

t = 2.00

480.578

2412.466

0.199

t = 2.50

486.638

2505.145

0.194

may still be achieved in the concentrations at 4.00, 6.00, 8.00, and 10.00 h (see Figure 6i,j,k,l). Here, the 6.00 h time point after dosing is the one showing a multimodal distribution with the highest standard deviation. Note that in the third column of Table 3 the highest standard deviation is 308.132, which corresponds to 6.00 h after the drug is administered. Therefore, simulated data is added for this time point. The Bayesian approach is again applied now considering two data points (i. e., 1.50 h and 6.00 h) and the population prior. The distributions of drug concentration resulting from using the new posterior distribution of PK parameters are presented in Figure 7. Figure 8 compares the probability distributions of drug concentration obtained by using (i) the population prior and (ii) the probability distributions after having samples at 1.50 h and 6.00 h. As can be noticed, the variability can be remarkably reduced and the PK parameters can be well characterized by having only two samples. Thus, the proposed sampling schedule for new patients is to have blood samples at 1.50 h and 6.00 h

after dosing. This also demonstrates the capabilities of the Bayesian approach for discriminating a patient given the intersubject PK variability information by means of the population prior P . In the following section, this capability is exploited for determining the PK parameters for new patients. Parameter Estimation for New Patients. In order to demonstrate the capabilities of the proposed approach, subjects P01, P06, and P10 are assumed to be new patients. We would like to emphasize that the population prior for an assumed new patient case is built by using the PK parameters of all other patients included in the clinical study except those of that patient to properly assess the advantages of the proposed approach. The new patient PK parameters are estimated again by applying the Bayesian approach considering (i) the previously proposed schedule for blood samples (i.e., 1.50 h and 6.00 h) and (ii) the corresponding population prior as explained in the Stage II section. Notice that the selected sampling points are indeed included in the available clinical data; otherwise, the closest ones to each of the selected experimental points would have been used for the purposes of evaluating the capabilities of this approach. The posterior marginal probability distributions of the PK parameters for the new patient P01 are presented in Figure 9. Similar results were obtained for the other assumed new patients. One can see how the variability has been reduced in the predicted response (i.e., drug concentration) from the population prior in 5123



ARTICLE

Table 3. Summary for the Predictive Concentration Distributions at Different Times after Dosing for Different Candidates for First Sample Point candidates for first sampling point sample point = 1.50 h time after dosing (h)

sample point = 2.00 h

mean (ng/mL)

SD (ng/mL)

CV

0.75

928.958

187.371

1.00

1238.722

192.429

1.25 1.50

1487.864 1686.620

1.75 2.00

sample point = 2.50 h

mean (ng/mL)

SD (ng/mL)

CV

mean (ng/mL)

SD (ng/mL)

CV

0.202

933.339

221.923

0.155

1241.353

224.945

0.238

978.994

271.599

0.277

0.181

1289.728

283.972

188.036 179.024

0.126 0.106

1487.499 1682.878

0.220

212.244 191.946

0.143 0.114

1533.055 1722.467

271.641 246.174

0.177 0.143

1843.322

169.753

0.092

1964.823

163.526

0.083

1836.289

170.595

0.093

1868.268

215.944

0.116

1954.823

153.477

0.079

1978.503

187.266

2.50

2124.031

166.555

0.095

0.078

2109.404

145.454

0.069

2116.569

154.027

3.00

2199.467

0.073

187.620

0.085

2181.931

170.176

0.078

2174.455

165.191

3.50

0.076

2216.047

216.177

0.098

2196.997

207.179

0.094

2177.512

201.543

0.093

4.00

2191.598

244.246

0.111

2172.045

243.017

0.112

2143.215

241.455

0.113

6.00 8.00

1893.047 1515.296

308.132 306.618

0.163 0.202

1876.089 1502.098

322.125 322.223

0.172 0.215

1830.308 1457.226

336.982 343.672

0.184 0.236

10.00

1177.634

277.251

0.235

1167.164

290.310

0.249

1129.076

312.195

0.277

12.00

904.165

241.529

0.267

895.519

251.342

0.281

865.092

270.625

0.313

14.00

691.140

207.732

0.301

683.800

214.745

0.314

660.167

230.299

0.349

18.00

403.909

152.301

0.377

398.490

155.723

0.391

384.651

164.281

0.427

24.00

183.866

94.194

0.512

180.555

95.436

0.529

174.322

98.104

0.563

30.00

86.206

56.821

0.659

84.306

57.377

0.681

81.368

57.917

0.712

195.514

0.330

201.604

0.348

223.274

0.366

average

Figure 6. Drug concentration (ng/mL) distribution for different time points by using the population prior and simulated concentration data for 1.50 h after dosing.

Figure 10. This figure shows the predicted concentration distributions at different time points for both the population prior and the individualized distribution obtained for the

patient P01. Additionally, the 95% highest probability density (HPD) bands of blood concentration for the population prior and the individualized distributions are compared in Figure 11. 5124



ARTICLE

Figure 7. Predicted drug concentration (ng/mL) distribution for different time points by using the population prior and simulated concentration data for 1.50 h and 6.00 h after dosing.

Figure 8. Comparison of predicted drug concentration (ng/mL) distributions for different time points (solid line, posterior distribution using two sample points (1.50 h and 6.00 h); dashed line, population prior).

Note how the variability of individualized posterior distributions for new patients (by applying the Bayesian approach using a new patient’s data from just two time points) is reduced in comparison to the population variability. Dosage Regimen Optimization. Once the joint posterior probability distribution of PK parameters is obtained for the new patient, the dosage regimen can be determined so that drug concentration is kept within the therapeutic window. Here, it is assumed that the therapeutic window is given and fixed. The proposed therapeutic window for Gabapentin to be used as a seizure control agent is 2-10 μg/mL.36,37

Figure 9. Marginal posterior distribution of PK parameters for subject P01 (solid line, patient distribution; dashed line, population prior).

According to Parke-Davis,38 the effective dosage of Gabapentin for seizure treatment is 900-1800 mg/day given in three equal doses. A starting dose of 300 mg, three times a day is recommended. If it is necessary, the dose may be increased up to 1800 mg/day. Table 4 shows the confidence levels obtained if these nominal dosage regimens are followed by the patients under study. For example, notice that the starting dose of 300 mg, three times a day renders a 49.7% confidence that the concentration of the drug will be within the therapeutic window during the 8 h interval (i.e., global approach) for patient P10. Furthermore, there is no new patient, except P10, for whom a concentration confidence band can be confined to the therapeutic window at 5125



ARTICLE

Figure 10. Predicted drug concentration (ng/mL) distribution for different time points by using the population prior and simulated data for 1.50 h and 6.00 h after dosing for subject P01 (solid line, patient’s distributions; dashed line, population prior).

any time during the dosing interval (i.e., local approach) by following any of the recommended dosages. Moreover, the confidence level obtained by administering 600 mg, three times a day to patient P10 is very poor, just 54%. Note that the recommended dosage, which is determined by the pharmaceutical companies adopting a “one dose fits all” approach, may not be the most appropriate regimen for every individual. Specifically, this may lead to poor Gabapentin therapy results for the patients studied since the concentration may fall below the efficacy level. Figure 12 shows a 95% drug concentration confidence band at the steady state for patient P06 by following the recommended dosages. Table 5 shows the optimized dosage regimens for the different patients which have been determined by employing the methodology described in the Stage III: Dosage Regimen Optimization section. For comparison purposes, this table also includes the dosage regimen that would be obtained if the total data set available for a patient is used. It is noteworthy that by just reducing the dosing interval to 6 h and administering the proposed 600 mg for patient P01 and P10, we can obtain a 95% confidence level that the concentration can be placed within the therapeutic window. However, a dosage regimen of 500 mg every 4 h is a more suitable therapy for patient P06; again by doing so, the concentration confidence band with a confidence level greater than or equal to 95% is confined to the therapeutic window (see Figure 13). These examples demonstrate how the recommended regimen (usually optimized for the population) may differ from the most adequate dosage regimen for an individual. Moreover, they emphasize the significance of a customized dosage regimen for providing a safer and more effective therapy. We would like to point out that one of the aims of this work is to exploit the information that is underlying the data derived from clinical trials. The idea is to use the population prior (P ) to capture this information. Thus, our methodology relies on the quality of the clinical trial as input. There are statistical methodologies (e.g., parallel group designs, cluster randomized designs, crossover designs) which have been proposed for the design and analysis of clinical trials and sample size determination. Such methodologies intend to support this process so that the studied subjects provide a reasonable representation of the

population. Adaptive and sequential approaches have also been developed which re-estimate the trial size as new data is being collected. Consequently, the sample size may change in the course of a clinical trial and not be a fixed number given a priori. Methodologies to determine the adequate trial size are beyond the scope of this paper; however, there is an extensive body of literature dealing with this issue.39 For the Gabapentin case study that we have presented in this work, the available data was limited to 19 patients (the clinical study executed to gather this data was not for FDA approval purposes). Nonetheless, it allowed us to show the benefits of the proposed Bayesian framework. Computational Aspects. The Bayesian analysis was performed using the R software40 and the MCMCpack package.41 It is noteworthy that R is free software under the terms of the GNU General Public License. The main advantage of this platform is that it provides a general purpose Metropolis sampling algorithm which facilitates the definition of customized prior probability distributions and likelihood functions. This has been especially important when the prior distribution is given by the prior population (P ). In stage I, the posterior PK parameters’ distribution of a patient included in the clinical trial has been obtained with a sample size of 3.0 105 parameter observations. It takes on average 222.4 CPU seconds to get a posterior distribution for such a case. On the other hand, in stage II, it takes on average 4828.4 CPU seconds to obtain the posterior distribution for a new patient (i.e., the population prior P is being used) having the same sample size. Convergence was already achieved after 3.0 105 iterations for both cases. The convergence was evaluated using cumulative plots for the average likelihood of the chains and trace plots or sample paths for the parameters.42 All computations have been carried out on an Intel i5 2.66 GHz computer.

’ FINAL REMARKS This work proposes a Bayesian approach for customized health care. Specifically, we propose a methodology for determining an individualized dosage regimen which becomes relevant for drugs whose PK varies widely among patients, adverse reactions are severe, or the therapeutic window is narrow. 5126



ARTICLE

Table 4. Confidence Levels Obtained by Administering the Standardized Dosage proposed Bayesian approach dosing interval

dose (mg)

local

global

patient P01 8h

300

n/da

48.19%

8h

600

n/d

83.80%

8h

300

n/d

38.60%

8h

600

n/d

85.34%

8h

300

n/d

49.65%

8h

600

54%

96.70%

patient P06

patient P10

a

There is no confidence level that places its corresponding concentration confidence band within the therapeutic window.

Figure 12. A 95% predicted drug concentration confidence band at steady state for patient P06 by following the standardized dosage regimen.

Figure 11. 95% HPD predicted blood concentration bands for new patients (the points represent the experimental data utilized to compute the individualized posterior distribution).

A Bayesian study relies on two elements: (i) the prior knowledge and (ii) experimental data. Here, we suggest exploiting the pharmacokinetic data collected during clinical trials to build a

population prior distribution (P ). We also described how, by using this population prior, the number of experiments required to characterize the pharmacokinetic parameters of a new individual can be minimized. We formally present the dosage regimen optimization problem and show how this problem can be simplified by fixing the dosing interval. This is a suitable approach for determining a dosage regimen whose dosing interval can easily be followed by the patient. The proposed approach has been applied to data obtained during a clinical study of Gabapentin pharmacokinetics. It has been demonstrated that the proposed Bayesian approach can characterize patients by using very few samples, two for the Gabapentin case, given the population prior. This case study has also demonstrated that the nominal dosage regimen recommended by pharmaceutical companies may not be the most appropriate therapy for an individual to follow due to the interpatient variability. It has been shown that for some patients the recommended therapies may provide poor results; namely, concentrations are not guaranteed to be confined within the therapeutic window. However, an optimized dosage regimen can be obtained by using the individualized distributions that result from applying the Bayesian approach presented herein. The optimized dosage regimen results in a drug concentration profile which is confined to the therapeutic window given a target confidence level, thus providing a safer and more effective therapy and very likely reducing the health care costs. Our future work is focused on the inclusion of covariates such as genomics, age, and gender in Bayesian modeling. Also, the 5127



ARTICLE

Table 5. Individualized Dosage Regimens for a Confidence Level of 95% feasible dose range (mg) proposed Bayesian approach dosing interval (h)

local

4.0

[271.9, 591.9]

6.0

total data set

global

local

global

[219.5, 611.3]

[277.2, 571.0]

[221.8, 597.1]

[573.5, 821.9]

[388.8, 860.2]

[635.3, 789.0]

[406.6, 835.9]

8.0

n/da

[661.0, 1066.6]

n/d

[725.0, 1029.7]

12.0

n/d

n/d

n/d

n/d

4.0

[346.0, 619.9]

[272.4, 648.3]

[484.3, 608.3]

[342.3, 643.1]

6.0

n/d

[533.2, 908.5]

n/d

[659.6, 913.3]

8.0

n/d

[985.5, 1108.8]

n/d

n/d

12.0

n/d

n/d

n/d

n/d

4.0

[278.0, 573.7]

[225.6, 593.6]

[216.3, 663.9]

[201.4, 685.7]

6.0

[581.2, 802.8]

[392.4, 842.4]

[411.0, 912.4]

[363.2, 952.1]

8.0

n/d

[641.6, 1046.8]

[752.9, 1096.8]

[618.8, 1150.4]

12.0

n/d

n/d

n/d

n/d

patient P01

patient P06

patient P10

a

There is no feasible dose range for a confidence level of 95%.

’ NOTATION Sets i = experimental observations j = patients k = subset of new patients Parameters

Figure 13. 95% predicted concentration confidence band at the steady state for an optimized dosage regimen for patient P06; τ = 4 h, Dose = 500 mg.

utilization of optimization approaches to obtain the individualized distributions are worthy to explore to reduce the considerable computational burden associated with the MCMC approach. In this work, we have assumed that the therapeutic window is fixed and the same for every individual. Individualized pharmacodynamic models are required to relax this assumption, which will be further investigated in our program.

’ AUTHOR INFORMATION Corresponding Author

*E-mail: [email protected].

’ ACKNOWLEDGMENT Support from the United States National Science Foundation (Grant NSF-CBET-0941302) is gratefully acknowledged. We would like to thank University of California, San Francisco for providing the data used in this work.

Dj = vector of experimental data for patient j Dose = dose amount Dosemin = minimum allowed value for the dose amount F = bioavailability ka = absorption rate ke = elimination rate to = dead time TWL = lower limit for the therapeutic window (i.e., efficacy level) TWU = upper limit for the therapeutic window (i.e., toxicity level) uj = vector of controlled variables for patient j in the general PK model V = volume of distribution zj = vector of covariates for patient j R = confidence level γ = an error parameter in the error model ζ = an error parameter in the error model ϑ = vector of unknown error parameters τ = dosing interval τmin = minimum allowed value for the dosing interval φ = vector of unknown PK parameters ψ = combined vector of PK and error parameters ω = an error parameter in the error model Ω = variance-covariance matrix of between subject variability Variables ^ t = predicted drug concentration at time t C xj = vector of state variables for patient j in the general PK model 5128


Industrial & Engineering Chemistry Research ^yj = vector of predicted values for the observations for patient j in the general PK model Functions

f = vector function relating the state variables xj to the predicted observations F = a generic multivariate probability distribution g = vector of differential equations that define a general PK model h = vector of algebraic equations that define a general PK model L(Dj|φ,ϑ) = likelihood of data Dj given a specific set of parameters N = an univariate/multivariate normal distribution p(φ,ϑ) = PK parameters prior distribution p(φ,ϑ|Dj) = PK parameters posterior distribution for patient j given the experimental data Dj P = the population prior Q φ ðt, τÞ = function to define concentrations for multidose models (dosage regimen optimization - local approach) φ ðtÞ = function to define concentrations during dosing interval Q τ for multidose models (dosage regimen optimization —global approach) ε = experimental error σ = vector of standard deviations for the experimental error

’ REFERENCES (1) Rehnquist, J. FDA’s Review Process for New Drug Applications: A Management Review; Technical Report OEI-01-01-00590, 2003. (2) Williams, P. J.; Ette, E. I. Pharmacometrics: Impacting Drug Development and Pharmacotherapy. In Pharmacometrics: The Science of Quantitative Pharmacology; Ette, E. I., Williams, P. J., Eds.; WileyInterscience: Hoboken, NJ, 2007; Chapter 1, pp 1-21. (3) Eichelbaum, M.; Ingelman-Sundberg, M.; Evans, W. E. Pharmacogenomics and individualized drug therapy. Ann. Rev. Med. 2006, 57, 119–137. (4) Krahenbuhl-Melcher, A.; Schlienger, R.; Lampert, M.; Haschke, M.; Drewe, J.; Krahenbuhl, S. Drug-related problems in hospitals - A review of the recent literature. Drug Safety 2007, 30, 379–407. (5) Steed, J. Drugged-out Seniors a Prescription for Disaster, Toronto Star, 2008. (6) Boissel, J. P.; Durieu, I.; Girard, P.; Nony, P.; Chauvin, F.; Haugh, M. Dose-Ranging Trials: Guidelines for Data Collection and Standardized Descriptions. Controlled Clinical Trials 1995, 16, 319–330. (7) Brodie, B. B. Physicochemical and Biochemical Aspects of Pharmacology. JAMA, J. Am. Med. Assoc. 1967, 202, 600–609. (8) Buxton, I. Pharmacokinetics and Pharmacodynamics: The Dynamics of Drug Absorption, Distribution, Action, and Elimination. In Goodman & Gilman’s The Pharmacological Basis of Therapeutics; Brunton, L. L., Lazo, J. S., Parker, K. L., Eds.; McGraw-Hill: New York, 2006; Chapter 1. (9) Wagner, J. G. Pharmacokinetics for the Pharmaceutical Scientist; Technomic Publishing Company: Lancaster, PA, 1993. (10) Blau, G.; Lasinski, M.; Orcun, S.; Hsu, S.-H.; Caruthers, J.; Delgass, N.; Venkatasubramanian, V. High fidelity mathematical model building with experimental data: A Bayesian approach. Comput. Chem. Eng. 2008, 32, 971–989. (11) Ette, E. I.; Williams, P. J.; Ahmad, A. Population Pharmacokinetic Estimation Methods. In Pharmacometrics: The Science of Quantitative Pharmacology; Ette, E. I., Williams, P. J., Eds.; Wiley-Interscience: Hoboken, NJ, 2007; Chapter 10, pp 265-285. (12) Ette, E. I.; Williams, P. J. Population Pharmacokinetics I: Background, Concepts, and Models. Ann. Pharmacother. 2004, 38, 1702–1706.

ARTICLE

(13) Reilly, P. M.; Blau, G. The use of statistical methods to build mathematical models of chemical reaction systems. Can. J. Chem. Eng. 1974, 52, 289–299. (14) Ruppert, D.; Cressie, N.; Carroll, R. J. A Transformation/ Weighting Model for Estimating Michaelis-Menten Parameters. Biometrics 1989, 45, 637–656. (15) Chan, P. L. S.; Weatherley, B.; McFadyen, L. A population pharmacokinetic meta-analysis of maraviroc in healthy volunteers and asymptomatic HIV-infected subjects. Br. J. Clin. Pharmacol. 2008, 65, 76–85. (16) Duffull, S. B.; Friberg, L. E.; Dansirikul, C. Bayesian Hierarchical Modeling with Markov Chain Monte Carlo Methods. In Pharmacometrics: The Science of Quantitative Pharmacology; Ette, E. I., Williams, P. J., Eds.; Wiley-Interscience: Hoboken, NJ, 2007; Chapter 5, pp 137164. (17) FDA, Guidance for Industry: Population Pharmacokinetics; U.S. Department of Health and Human Services: Washington, DC, 1999. (18) Jelliffe, R.; Schumitzky, A.; Guilder, M. V.; Liu, M.; Hu, L.; Maire, P.; Gomis, P.; Barbaut, X.; Tahani, B. Individualizing drug dosage regimens: roles of population pharmacokinetic and dynamic models, Bayesian fitting, and adaptive control. Ther. Drug Monit. 1993, 15, 380– 393. (19) Steimer, J.-L.; Mallet, A.; Golmard, J.-L.; Boisvieux, J.-F. Alternative Approaches to Estimation of Population Pharmacokinetic Parameters: Comparison with the Nonlinear Mixed-Effect Model. Drug Metab. Rev. 1984, 15, 265–292. (20) Sheiner, L.; Grasela, T. J. Introduction to mixed effect modeling: concepts, definitions, and justification. J. Pharmacokinetics and Biopharmaceutics 1991, 19, 11S–24S. (21) Mallet, A. A maximum likelihood estimation method for random coefficient regression models. Biometrika 1986, 73, 645–656. (22) Davidian, M.; Gallant, A. R. Smooth nonparametric maximum likelihood estimation for population pharmacokinetics, with application to Quinidine. J. Pharmacokinet. Biopharm. 1992, 20, 529–556. (23) Whiting, B.; Kelman, A. W.; Bryson, S. M.; Derkx, F. H.; Thomson, A. H.; Fotheringham, G. H.; Joel, S. E. Clinical pharmacokinetics: a comprehensive system for therapeutic drug monitoring and prescribing. Br. Med.l J. 1984, 288, 541–545. (24) Mehvar, R. Pharmacokinetic-based design and modification of dosage regimens. Am. J. Pharm. Educ. 1998, 62, 189–195. (25) Duffull, S. B.; Kirkpatrick, C.; Green, B.; Holford, N. Analysis of population pharmacokinetic data using NONMEM and WinBUGS. J. Biopharm. Stat. 2005, 15, 53–73. (26) Salinger, D. H.; McCune, J. S.; Ren, A. G.; Shen, D. D.; Slattery, J. T.; Phillips, B.; McDonald, G. B.; Vicini, P. Real-time Dose Adjustment of Cyclophosphamide in a Preparative Regimen for Hematopoietic Cell Transplant: A Bayesian Pharmacokinetic Approach. Clin. Cancer Res. 2006, 12, 4888–4898. (27) McCune, J. S.; Batchelder, A.; Guthrie, K. A.; Witherspoon1, R.; Appelbaum, F. R.; Phillips, B.; Vicini, P.; Salinger, D. H.; McDonald, G. B. Personalized Dosing of Cyclophosphamide in the Total Body Irradiation-Cyclophosphamide Conditioning Regimen: A Phase II Trial in Patients With Hematologic Malignancy. Clin. Pharm. Ther. 2009, 85, 615–622. (28) Bishop, C. M. Pattern Recognition and Machine Learning; Springer Science: New York, NY, 2006. (29) Barrett, J.; Mondick, J.; Narayan, M.; Vijayakumar, K.; Vijayakumar, S. Integration of modeling and simulation into hospital-based decision support systems guiding pediatric pharmacotherapy. BMC Med. Inf. Decis. Making 2008, 8, 6. (30) Salinger, D. H.; Vicini, P.; Blough, D. K.; O’Donnell, P. V.; Pawlikowski, M. A.; McCune, J. S. Development of a population pharmacokinetics-based sampling schedule to target daily intravenous busulfan for outpatient clinic administration. J. Clin. Pharmacol. 2010, 50, 1292-1300. (31) Hsu, S.-H.; Stamatis, S. D.; Caruthers, J. M.; Delgass, W. N.; Venkatasubramanian, V.; Blau, G. E.; Lasinski, M.; Orcun, S. Bayesian 5129



ARTICLE

Framework for Building Kinetic Models of Catalytic Systems. Ind. Eng. Chem. Res. 2009, 48, 4768–4790. (32) Metropolis, N.; Rosenbluth, A. W.; Rosenbluth, M. N.; Teller, A. H.; Teller, E. Equation of State Calculations by Fast Computing Machines. J. Chem. Phys. 1953, 21, 1087–1092. (33) Simonoff, J. S. Smoothing Methods in Statistics; Springer-Verlag: New York, 1996; Springer Series in Statistics. (34) Jeliffe, R.; Schumitzky, D.; Bayard, D.; Leary, R.; Thompson, A.; Boten, A.; Van-Guilder, M.; Bustad, A.; Neely, M. Pharmacokinetic tools to optimize control of drug PK/PD models for best patient care; Technical Report, 2008. (35) Urban, T.; Brown, C.; Castro, R.; Shah, N.; Mercer, R.; Huang, Y.; Brett, C.; Burchard, E.; Giacomini, K. Effects of Genetic Variation in the Novel Organic Cation Transporter, OCTN1, on the Renal Clearance of Gabapentin. Clin. Pharmacol. Ther. 2008, 83, 416–421. (36) ARUP-Laboratories, ARUP’s Laboratory Test Directory Gabapentin. http://www.aruplab.com/guides/ug/tests/0090057.jsp (accessed Jan 2011). (37) García, L. L.; Shihabi, Z. K.; Oles, K. Determination of gabapentin in serum by capillary electrophoresis. J. Chromatogr., Sect. B 1995, 669, 157–162. (38) Neurontin (gabapentin); Technical Report, Parke-Davis: Ann Arbor, MI, 2009. (39) Chow, S.-C.; Liu, J.-P. Design and Analysis of Clinical Trials: Concepts and Methodologies, 2nd ed.; Wiley Series in Probability and Statistics; Wiley-Interscience: Hoboken, NJ, 2003. (40) R Development Core Team. R: A Language and Environment for Statistical Computing; R Foundation for Statistical Computing: Vienna, Austria, 2008, ISBN 3-900051-07-0. (41) Martin, A. D.; Quinn, K. M.; Park, J. H. Markov chain Monte Carlo (MCMC) Package; Washington University in St. Louis: St. Louis, MO; University of California at Berkeley: Berkeley, CA, 2010. (42) Cosma, I. A.; Evers, L. Markov Chains and Monte Carlo Methods; African Institute for Mathematical Sciences: Cape Town, South Africa; University of Cambridge: Cambridge, U.K.; University of Glasgow: Glasgow, Scotland, 2010.

5130


Pharmacokinetic Based Design of Individualized ... - ACS Publications

Recommend Documents