Influence of Age and Gender on the Pharmacokinetics and Pharmacodynamics of Remifentanil: I. Model Development

After obtaining Institutional Review Board approval and informed consent, 65 volunteers classified as American Society of Anesthesiology (ASA) physical status 1 and 2 were enrolled in the trial. Our goal in volunteer recruitment was to have nearly equal representation in three age groups: young (20–40 y), middle-aged (40–65 y), and elderly (> 65 y). The ages ranged from 20–85 y. Thirty-eight men and 27 women were studied. Remifentanil blood concentration and EEG effect data were collected in four separate phases. In the first phase, 10 healthy men classified as ASA physical status 1 whose ages ranged from 22 to 38 y received a constant-rate infusion of remifentanil at doses ranging from 1 to 8 micro gram [centered dot] kg sup -1 [centered dot] min sup -1 for 20 min. [5]In the next phase, ten more healthy men (ages 23 to 39 y) were enrolled in a cross-over study to investigate the comparative pharmacokinetics and pharmacodynamics of alfentanil and remifentanil. [6]In the remifentanil arm of the trial, the volunteers received an infusion of 2 or 3 micro gram [centered dot] kg sup -1 [centered dot] min sup -1 for 10 min. [6]In a third phase, ten healthy women volunteers (ages 20 to 33 y) were studied according to the same study protocol. In a fourth phase, 35 healthy men and women (ages 40 to 85 y) were enrolled. All volunteers in this phase received a remifentanil infusion at 3 micro gram [centered dot] kg sup -1 [centered dot] min sup -1 until there was no further change in the EEG, as measured by the spectral edge frequency (infusion duration, 4–15 min). We studied 15 more persons using the same protocol (11 men and 4 women between the ages of 41 and 84 y) from whom we collected EEG but not concentration data.

Remifentanil was administered through an 18- or 16-G intravenous catheter inserted into a forearm vein. A 20-G arterial catheter was inserted into the radial artery for sample collection. During all phases, frequent 3-ml arterial blood samples were taken at preset intervals until 4 h after the end of the infusion. The rapid metabolism of remifentanil by nonspecific blood and tissue esterases required that all esterases be inactivated immediately after the arterial samples were taken. This was achieved by adding acetonitrile to the sample to denature the plasma proteins. Subsequently, methylene chloride was added to the mixture to extract the remifentanil into the organic phase, which was then separated and stored at -70 degrees Celsius until assay. The remifentanil analyses were performed using a high-resolution, gas chromatographic-mass spectrometry assay with a quantitation limit of 0.1 ng/ml and a paired aliquot coefficient of variation of less than 15% for concentrations greater than 0.1 ng/ml. [7,8] 

An initial exploratory two-stage pharmacokinetic analysis was performed to determine whether a two- or three-compartment mammillary model best described each individual's data. To develop the final population pharmacokinetic model, we used the approach of Mandema and colleagues. [9]Two- and three-compartment population pharmacokinetic models without covariates were estimated with NONMEM* using the “first order conditional estimation” method and “eta-epsilon interaction” to reduce the influence of model misspecification. The interindividual error on each of the model parameters (V¹, V², V³, Cl¹, Cl sub 2, and Cl³) was modeled using a log-normal variance model:Equation 1where Pⁱis the value of the parameter in the individual, theta^TVis the typical value of the parameter in the population, and eta is a random variable with mean zero and variance omega sup 2. Technically, omega is the standard deviation of eta in the log domain. However, when omega is small, it is approximately the coefficient of variation (%CV) of P. To simplify the interpretation of our results, we report the variability of the pharmacokinetic parameters, P, as “%CV,” understanding that this is an approximate interpretation of the modeled variability. A “constant c.v.” model was used for the residual intraindividual error.

Empirical Bayesian [10,11]estimates of the individual pharmacokinetic parameters were obtained based on the typical values and variances of the estimated volumes and clearances. Unit disposition functions were calculated from these Bayesian estimates and plotted between the times of the first and last observations. To confirm that the population parameters provided a good description of the individual unit disposition functions, we superimposed the unit disposition function for the typical person on this graph. The median of the weighted residuals (calculated as the difference between the individual and the population unit disposition functions at the time of each sample) was computed as a measure of the bias over the first 10 min.

A generalized additive model (GAM) was used to permit identification of linear and nonlinear relations between the Bayesian estimates of individual volumes and clearances and the volunteer covariates. We provided the GAM function (as implemented in S-PLUS**) with a list of covariates available for analysis. The covariates analyzed were age, gender, weight, height, body surface area, and lean body mass (LBM). Body surface area was calculated from the weight (measured in kilograms) and height (measured in centimeters) using the formula of DuBois and DuBois [12]:Equation 2. Lean body mass was calculated from the gender, weight (in kilograms), and height (in centimeters)***:Equation 3, Equation 4.

For each of the covariates in this list, we provided the GAM function with an ordered regimen of possible candidate forms describing how each covariate could enter the model. We allowed three forms for each covariate. Each covariate could (1) be omitted from the model;(2) enter the model linearly; or (3) enter the model nonlinearly. We chose a natural cubic spline with knots at the 33rd and 66th percentile of the data as a flexible function to search for potential nonlinear relations. A nonlinear form was not included for gender. The GAM function performed a stepwise search to find the significant covariates for each of the volumes and clearances and to find the best form (linear or nonlinear) of each important covariate. Not all possible combinations of covariates and forms were examined. The model that resulted in the biggest decrease in the Akaike Information Criterion was returned by the GAM function as the best model. Ninety-five models were examined in this step. The GAM function is described in more detail by Hastie. [13] 

The structural model incorporating covariates identified by the GAM analysis was then evaluated with NONMEM to develop the final parameter estimates. The final model was examined for parsimony by excluding individual covariates and demonstrating a statistically significant increase in the NONMEM objective function, and by analysis of the standard errors of the estimated parameters.

As previously defined, we used the weighted residual to describe the quality of the prediction in each individual, calculated as (measured-predicted)/predicted. [14]We examined the performance of the complex pharmacokinetic model in the worst case for young, middle-aged, and elderly subjects. These persons were selected based on the worst mean absolute weighted residual in each age group. We also examined the quality of the prediction of the pharmacokinetic models for the population by calculating the median absolute weighted residual (MDAWR) and the mean of the individual mean absolute weighted residuals (MAWR) as described by Kataria and associates, [14]and by examining the pattern of the weighted residuals, expressed as measured/predicted value. [15] 

In anticipation of developing dosing guidelines based on weight, we examined the performance of a weight proportional model and a lean body mass proportional model, in which all of the volumes and clearances are assumed to increase in proportion to the individual's weight or lean body mass, respectively.

The EEG was digitized and the spectral edge calculated as previously reported. [16]We performed an initial two-stage pharmacodynamic analysis, calculating the effect site concentration using the individual two- or three-compartment pharmacokinetic parameters obtained during the two-stage pharmacokinetic analysis. To develop the final population pharmacodynamic model, we again used the approach of Mandema and colleagues. [9]NONMEM's first-order conditional estimation method was used to analyze the pharmacodynamic data (spectral edge frequency) with a sigmoid E^maxmodel. Equation 5where E⁰is the measured baseline effect, E^maxis the maximum measured effect, Ce is the apparent effect site concentration, EC⁵⁰is the concentration that causes 50% of the maximum effect, and gamma describes the steepness of the concentration response relation. The apparent effect site concentration was calculated for each individual based on the individual dosing information, previous Bayesian estimates of the individual pharmacokinetic parameters (which entered into the pharmacodynamic model as known parameters or “fixed effects”), and the equilibration rate constant, k^e0(estimated by NONMEM). We used a log-normal variance model to describe the residual interindividual error on each of the model parameters (k^e0, E⁰, E^max, EC⁵⁰, gamma) and an additive model to describe the residual intraindividual error (in contrast to the constant c.v. model used for the pharmacokinetic analysis). The relation between the Bayesian estimates of the individual pharmacodynamic parameters and the volunteer covariates was also examined using a GAM, as described previously. Seventy models were evaluated in this step.

As described for the pharmacokinetic analysis, the pharmacodynamic model incorporating covariates identified by the GAM analysis was evaluated with NONMEM to develop the final pharmacodynamic model. The final model was examined for parsimony by excluding individual covariates and demonstrating a statistically significant increase in the NONMEM objective function, and by analyzing the standard errors of the estimated parameters.

In contrast to the kinetic analysis, we used the unweighted residual (observed-predicted) to describe the quality of the prediction in each individual. We examined the performance of the simple and complex pharmacodynamic model in the worst case for young, middle-aged, and elderly persons. We also examined the quality of the prediction of the pharmacodynamic models for the population by calculating the median absolute residual and the mean of the individual mean absolute residuals (MAR). To provide some sense of the magnitude of these values, we also reported them as a percentage of the pharmacodynamic range (E⁰-E^max) of the respective model.

We used the dose-versus-response data from the 15 additional participants to compare prospectively the performance of the simple (without covariates) and complex (with covariates) pharmacokinetic and pharmacodynamic models. We examined the predictive performance of the simple and complex pharmacokinetic and pharmacodynamic models in the the worst case for the middle-aged and elderly subjects. We examined the quality of the predictions of the combined model by calculating the median absolute prediction error and the mean of the individual mean absolute prediction errors (MAPE) as measures of the precision of the models. We also examined the median prediction error and the mean of the individual mean prediction errors (MPE) as measures of the bias of the models, as detailed by Varvel and coworkers. [17]However, we defined prediction error as the difference between the observed EEG effect and the predicted EEG effect, in keeping with our use of an additive variance model for the residual intrasubject error. These values were reported as a percentage of the pharmacodynamic range (E⁰-E^max) of the respective model.

All 65 participants enrolled completed the study. Table 1summarizes their distribution by age group and gender, and Figure 1shows the individual values by age, gender, weight, height, LBM, and body surface area.

The initial two-stage pharmacokinetic analysis found that the concentration-time profiles for approximately one half of the individuals were best described using a two-compartment model, and the other half were best described using a three-compartment model. The plots of V¹, Cl sub 1, and Vd^as(calculated as V¹+ V²for the two-compartment models, and as V¹+ V²+ V³for the three-compartment models) showed decreasing trends with increasing age.

The population model derived using NONMEM found that the pharmacokinetics of remifentanil were best described by a three-compartment model (P << 0.01). The parameter estimates and performance of the three-compartment model without covariates are shown in Table 2, together with the results of the weight proportional and the LBM proportional models. The LBM proportional model resulted in a greater improvement in the objective function than the weight proportional model, although this was associated with only a small improvement in the MDAWR and (MAWR) with bar. The main improvement resulting from the LBM proportional model was a 10% decrease in the %CV of V¹, together with small decreases in the %CV of V²and Cl¹.

(Figure 2) shows the unit disposition functions obtained with the simple three-compartment model. The unit disposition functions were calculated from the time of the first to the time of the last detectable concentration. Although sampling was performed for 4 h beyond the end of infusion, the remifentanil concentrations in many individuals were less than the limit of detection before the third exponential phase was apparent. The simple model (without covariates) described the early clinically important component of the unit disposition functions well, with a small negative bias (-5.75%) for the first 10 min.

The GAM method was used to identify potentially significant covariates. The best models included age for all volumes and clearances. In addition, LBM was a covariate of V¹, V², and Cl¹in the best models. These were the same parameters for which inclusion of LBM as a covariate decreased the %CV in the LBM-adjusted model described before. No significant nonlinear relations were detected. Figure 3shows the individual Bayesian estimates of V¹and Cl¹as a function of age.

An initial model with age on all parameters and LBM on V¹, V², and Cl¹was implemented in NONMEM. The covariates were modeled around their median value, such that the intercept term reflects the typical value for a person of median age and LBM. This initial model was further refined by deleting certain parameters according to the criteria outlined by Mandema and colleagues. [9]Only the linear relations between age and V³did not remain in the final model (Table 3). From the age of 20 to 85 yr, we found that typical values for V¹and Cl¹decreased by approximately 25% and 33%, respectively.

The more complex model accounted for a 9% to 16% decrease in the %CV of the structural model parameters. In addition, the standard errors on the variance of eta for all parameters except V³were decreased by 15%(V¹) to 50%(Cl¹). This more complex pharmacokinetic model showed a substantial improvement in -2 log likelihood (213.5; P < 0.01) and was associated with an improvement in the MDAWR and MAWR with bar of approximately 5% and 4%, respectively. The improvement was consistent across all age groups, as assessed by the worst mean absolute weighted residuals. Figure 4shows the weighted residuals for the simple model, the LBM proportional model, and the complex pharmacokinetic model. The improvement with the addition of both LBM and age as significant covariates is seen in the lower panel.

Four participants were excluded from the pharmacodynamic analysis because of very noisy or no EEG data. The initial two-stage pharmacodynamic analysis detected an age-related decrease in EC⁵⁰and k^e0. Figure 5shows the individual Bayesian estimates of the concentration-response relations, together with the typical response of the population model without covariates. This figure suggests that there is considerable interindividual variability in E⁰, E^max, EC⁵⁰, and gamma. Table 4shows the parameters and performance of the simple pharmacodynamic model. The median absolute residual for the simple model was 2.11 Hz (15% of the dynamic range).

The GAM method was used to identify potentially significant covariates. The best models included an age-related linear decrease in k sub e0 and EC⁵⁰. No significant nonlinear relations were detected. The decreases in k^e0and EC⁵⁰with age were confirmed by subsequent analysis with NONMEM. Table 5shows the final age-adjusted pharmacokinetic model. The parameters were modeled about the median age, such that the intercept reflects the typical value for a 40 y old. Figure 6illustrates the doubling of t¹/2 k^e0and the halving of EC sub 50 that occur from ages 20 to 85 y in the typical individual. The complex pharmacodynamic model showed a significant improvement in -2 log likelihood (27; P < 0.01) associated with only a small improvement in the median absolute residual and MAR with bar. This model resulted in an 8% decrease in the %CV of the k^e0and a 9% decrease in the %CV of the EC sub 50. As with the pharmacokinetic analysis, the standard errors on the estimates of the variance of eta were decreased for those parameters with significant covariates in the complex model. The standard errors of k^e0and EC⁵⁰were decreased by 17% and 11%, respectively. The worst median absolute variance in the young and middle-aged groups also improved. Some improvement in the pattern of the residuals was noted after infusion (Figure 7).

The 15 participants (11 men, 4 women) in whom the simple and complex pharmacokinetic and pharmacodynamic models were tested prospectively were ages 41 to 84 y and had LBMs of 41 to 74 kg. The median absolute prediction error of the simple model was 23% of the pharmacodynamic (E⁰-E^max) range. The median absolute prediction error of the complex model was 19% of the pharmacodynamic range. The improvement in median absolute prediction error in this prospective analysis suggests that the inclusion of age and LBM in the complex model will allow more accurate model performance. However, this improvement came at the cost of an increase in bias: The median prediction error was 0.5% with the simple model and 7% with the complex model.

Unique features of remifentanil are its rapid clearance and rapid k^e0, resulting in a rapid onset and offset of drug effect. It is tempting to speculate that these characteristics will make remifentanil an easy drug to titrate, and that clinicians will not need to consider patient covariates when choosing a dosing regimen. However, the rapid onset of drug effect may be accompanied by rapid onset of adverse events such as apnea and muscle rigidity. The rapid offset of drug effect can cause severe patient pain when the anesthesiologist is ill equipped to deal the problem, such as when the patient is being transported to the recovery room. In a companion manuscript, [18]we explore the dosing implications of the influence or age, lean body mass and intersubject variability on remifentanil dosing guidelines.

Previous studies have reported conflicting findings concerning the influence of age and gender on the pharmacokinetics of opioids. Helmers and associates [19]compared alfentanil kinetics between young and elderly patients after a single bolus administration. They found that the plasma clearance of alfentanil is less rapid and the elimination half-life is longer in elderly patients than in young ones. In contrast, Scott and Stanski [2]found no age-related changes in the pharmacokinetic parameters of fentanyl or alfentanil. Hudson and colleagues [20]found that patient age was positively correlated with both alfentanil Vd^asand elimination half-life. Lemmens and associates [1]considered that the effect of gender could be a possible factor to explain the differences among various studies. These authors studied the effects of age on the pharmacokinetics of alfentanil in female and male patients undergoing lower abdominal surgery with nitrous oxide-alfentanil anesthesia. They found a significant negative correlation between plasma alfentanil clearance and age in women but not in men and concluded that the effects of age on the pharmacokinetics of alfentanil depend on gender. Matteo and coworkers [21]found that the initial volume of sufentanil distribution was significantly smaller in elderly patients, but that the elimination half-lives, plasma clearances, and total volumes of distribution were similar for elderly and younger participants. These authors believed that age-related differences in the action of sufentanil could not be explained by the observed differences in the initial volume of distribution and that alterations in pharmacodynamics appear to be more important in the prolonged opioid effect seen in the elderly. Lehmann and colleagues [22]found no significant correlation among sufentanil's half-life, the volume of distribution, or clearance with age or body weight. In a more recent study, Helmers and colleagues [23]reported that the sufentanil plasma concentration-time curves for young adult and elderly participants could be superimposed. The elimination half-life, clearance, and volume of distribution did not differ between both groups. They concluded that the age-related differences in the action of sufentanil cannot be explained by its pharmacokinetic properties.

Conflicting pharmacodynamic findings for the opioids have also been reported. Scott and Stanski [2]found a 50% decrease in the dose requirement of fentanyl and alfentanil with increasing age. They found that brain sensitivity, as determined by EEG changes, increased significantly with age. Lemmens and coworkers [3,4]investigated the effects of age on the pharmacodynamics of alfentanil, as a supplement to nitrous oxide anesthesia, in 14 younger and 14 older women undergoing curative surgery for primary breast cancer. They could not demonstrate an effect of age on the plasma concentration-effect relations for response to intubation, skin incision, and the probability of needing naloxone at the end of surgery. Ausems and colleagues [24]examined the response of 37 patients between 23 and 54 y to perioperative stimuli using alfentanil and 66% nitrous oxide in oxygen. For the single events of intubation, skin incision, skin closure, and postoperative ventilation, each patient's clinical state was categorized as “responsive” or “nonresponsive” based on predefined criteria. The data from all patients were pooled and analyzed with logistic regression to define alfentanil concentrations required along with 66% nitrous oxide to provide a satisfactory state of general anesthesia. They also recorded multiple observations during surgery to determine the quantal response to stimulation in each patient. No relation between pharmacodynamic parameters and age were reported.

We evaluated several different pharmacokinetic models for the following reasons. The two-stage pharmacokinetic analysis identified both two- and three-compartment models, which made it difficult to average the parameters to obtain a two-stage population parameter set. In addition, one half of the individuals did not have the third-compartment pharmacokinetic parameter estimates necessary for the covariate analysis. The mixed-effects modeling approach resolved these difficulties. Bayesian estimation of individual pharmacokinetic and pharmacodynamic parameters resulted in shrinkage of the parameter estimates toward the population mean characteristic of population analysis while still providing the individual estimates necessary for the covariate analysis.

The weight-proportional models may offer advantages to the clinician because they enable dosing guidelines to be derived on a per-kilogram basis. However, both simple and weight-adjusted models must be evaluated before arbitrarily selecting per-kilogram dosing as the preferred dosing model. Lean body mass entered the final complex pharmacokinetic model in a time-dependent and nonlinear manner, and its influence was less important than age. The complex models enable us to develop dosing guidelines for clinicians that will account for the estimates of interindividual variability in pharmacokinetic and pharmacodynamic parameters. The complex models are probably the best choice for computer-controlled drug delivery, although they are also appropriate for developing dosing nomograms and specific guidelines for use in special populations.

Our study involved more participants than did the opioid studies we cited. Because many covariates were correlated (such as gender, height, weight, LBM, and body surface area), and the data set was large, an efficient search strategy was needed to relate the participant covariates to pharmacokinetic and pharmacodynamic parameters. The GAM procedure provided a method to efficiently identify significant covariates. It was also able to build a model in the presence of correlated covariates, which presented a problem for more graphical methods of covariate analysis. [25]For instance, in each age group, nearly all the men had a greater LBM than did the women (Figure 1). Simply plotting the individual estimates of V¹against gender revealed that V¹was significantly less in the women. Although gender enters into the calculation of LBM, the GAM procedure could detect that LBM was a more important covariate than gender. Even with this large data set, the GAM procedure was very fast. The entire GAM analysis for the pharmacokinetic and pharmacodynamic data was complete in approximately 2 min. In contrast, evaluation of the final models with NONMEM required several months of continuous processing on a 90-MHz Intel Pentium computer processor running Windows NT (Microsoft Corp., Redmond, WA) and an optimized 32-bit version of NONMEM (Intel, Santa Clara, CA). An additional benefit of the GAM procedure was that it provided good initial estimates for the NONMEM analysis. However, we do not claim that by using the GAM procedure and the subsequent testing of the models with NONMEM we have necessarily found the optimal model. There may be other models that we have not tested that describe the observations better than the models we described. Testing every model with NONMEM would require years of effort with the available resources, and the likely improvement (assuming better models exist) would be small given the good performance of the models reported. Therefore, we restricted our exploration to those models that were identified as promising by the GAM analysis.

In a smaller study of 24 patients, Westmoreland and coworkers [26]concluded that age did not influence the total clearance of remifentanil. In contrast, we found that age is an important covariate in the pharmacokinetics of remifentanil. We also found a reduced central volume of distribution in the elderly, which will result in higher initial blood concentrations after a given bolus dose than in younger persons.

Although the canonical univariate parameter for opioids reported by Gregg and colleagues [27]and Gambus and associates [28]is a more consistent measure of opioid drug effect than the spectral edge frequency (SE⁹⁵), we chose to use the latter for consistency with previous work investigating the influence of age on opioid-induced EEG changes. Based on the spectral edge as a measure of opioid drug effect, we reported an age-related increase in the half-life k^e0for remifentanil. This slower equilibration between blood and brain concentrations in the elderly will offset their higher initial blood concentrations. Consequently, peak effect site concentrations after a bolus dose will be similar in the elderly and in the young, although the time of the peak will be relatively delayed. Thus, based on the EEG model, there are no pharmacokinetic grounds for recommending reduced bolus doses in the elderly, and it is for pharmacodynamic reasons (the 50% reduction in EC⁵⁰in the elderly) that we recommend halving bolus doses in the elderly.

The maintenance infusion rate is calculated as the product of the target concentration and the clearance of the drug. Again, based on the EEG findings of a 50% decrease in EC⁵⁰in elderly participants, we can anticipate that there will be at least a 50% decrease in maintenance infusion rates in elderly persons. However, the clearance in the elderly is about 66% of clearance in younger persons. This further decreases the infusion rate required to maintain any given drug effect. The net effect of the change in remifentanil potency and clearance in elderly patients is that the infusion rate to maintain any level of drug effect is about one third the rate in a younger person (66% x 50%= 33%).

We did not use the same weighted measures of performance to evaluate the pharmacokinetic and pharmacodynamic data. The pharmacodynamic data spanned a range from approximately 5 to 20 Hz. In contrast to the pharmacokinetic data (Figure 2) this range is less than one order of magnitude. Within this range the variance of the residual intraindividual error was well described by an unscaled error model. Because we used an unscaled error model for the objective function, we believe it more appropriate to use an unscaled measure for the residual error. However, if we had used the same weighted measures of performance, we would have found that the MDAWR and the MAWR of the pharmacodynamic data were only slightly greater than for the pharmacokinetic data. Specifically, the MDAWR and MAWR were 15.3% and 22.8%, respectively, for the complex pharmacokinetic model (Table 3) and were 17.2% and 24.7%, respectively, for the complex pharmacodynamic model.

Our finding of an increased sensitivity in the elderly using the EEG model is consistent with the earlier work of Scott and Stanski [2]for fentanyl and alfentanil. Although the EEG is not, in and of itself, a clinically useful end point, it offers several experimental advantages over clinical measures of drug effect. The EEG is a nearly instantaneous measure of drug effect, a continuous measure of drug effect, and an objective and reproducible measure of drug effect. These benefits of using the EEG would be irrelevant if the EEG told us nothing about the time course of clinically important opioid drug effects or about the relative potency of the opioids. The EEG predicts a more rapid and evanescent effect from a bolus of remifentanil and alfentanil compared with fentanyl or sufentanil. This correlates with the observed time course of effect of these drugs in the operating room for such effects as ventilatory depression, sedation, and analgesia. Shafer and Varvel [29]reviewed the pharmacodynamics of fentanyl, alfentanil, and sufentanil and observed that the EEG predicted relative potencies among these drugs that were entirely consistent with the predicted EC⁵⁰or therapeutic ranges for clinical measures of drug effect. Hill and associates [30]examined analgesic response, side effects, and evoked potentials during steady-state opioid infusions and found a robust relation between plasma drug concentration and analgesic, ventilatory, and subjective-effect magnitudes for each opioid in the study. We are aware of no data that suggest that the relative potency for opioids as measured by the EEG differs from the relative potency based on clinical measures of drug effect. The opioid studies that did not observe changes in brain sensitivity with age, those of Lemmens and coworkers [3]and Ausems and colleagues, [24]used low-resolution pharmacodynamic end points and examined fewer participants. Consequently, their studies may have lacked the resolution and the statistical power to detect the changes we report. However, the EEG is not clinical anesthesia, and it is possible that changes in potency reflected in the EEG may not reflect changes in potency for some clinical end points.

Although statistically significant, the addition of age and LBM as covariates resulted in a relatively small improvement in the prediction of EEG effect when applied prospectively. This may have been due in part to the study design. Relatively high infusion rates were chosen to achieve maximum EEG effect in all participants for the purpose of developing the full pharmacokinetic and pharmacodynamic model. Thus the EEG rapidly changed from baseline to high amplitude, low-frequency delta waves in all participants. As neither E⁰nor E^maxwere related to subject covariates, it is not surprising that the simple and complex pharmacodynamic models performed similarly during the infusion. However, even with this study design, covariate-related differences in pharmacokinetics and pharmacodynamics should allow a more accurate prediction of the EEG effect during the postinfusion period. In the group used to derive the models (Figure 7) and the group used for the prospective test, the pattern of the residuals for the pharmacodynamic data improved during the postinfusion period, although this was relatively small. Thus, even when statistically significant relations between model parameters and participant covariates are detected, further analysis is required to evaluate the importance of these relations in the presence of known interindividual variability in model parameters.

We used NONMEM and GAM to develop simple and complex population pharmacokinetic and pharmacodynamic models, and we assessed these models prospectively. We found (1) an effect of age on the pharmacokinetic and pharmacodynamic parameters for remifentanil;(2) an effect of LBM on the pharmacokinetic parameters; and (3) no effect of gender on any pharmacokinetic or pharmacodynamic parameter. In our companion article, we will use the estimates of interindividual variability in pharmacokinetic and pharmacodynamic parameters provided by NONMEM to perform computer simulations that examine the influence of age, LBM, and intersubject variability on remifentanil dosing guidelines.

*Beal SL, Sheiner LB: NONMEM User's Guide. San Francisco, University of California San Francisco, 1979. The pharmacokinetic and pharmacodynamic models for the NONMEM and GAM analyses are available via the WWW at URL: http://pkpd.icon.palo-alto.med.va.-gov.

**S-PLUS for Windows, Version 3.2, copyright 1988, 1994 Mathsoft, Inc.

***James WPT: Research on obesity. London, Her Majesty's Stationary Office, 1976.