Pharmacokinetic-Pharmacodynamic Modeling of the Respiratory Depressant Effect of Alfentanil 

After we obtained institutional review board approval and written informed consent, we studied 14 men classified as American Society of Anesthesiologists physical status I or II who were undergoing major urologic surgery. Their median age was 42 yr (range, 19–72 yr) and their median weight was 82.5 kg (range, 68–108 kg).

The unpremedicated patients were studied before anesthesia was induced. After arrival in the induction room, standard monitoring (noninvasive blood pressure, electrocardiography, and pulse oximetry) devices were placed. Two intravenous (both forearms) and one arterial cannula (radial artery of the nondominant hand) were inserted. An intravenous infusion of alfentanil (2.3 [micro sign]g [middle dot] kg^-1[middle dot] min^-1) was started while the patients were breathing oxygen-enriched air (FI^O(2)= 0.5) over a tightly fitting continuous positive airway pressure mask. This enabled us to obtain respiratory rate, approximate minute ventilation, and end-expiratory partial pressure of carbon dioxide (PE^CO(2)) using the standard monitors of an anesthesia workstation (Cicero; Draeger, Luebeck, Germany). The infusion was discontinued when a cumulative dose of 70 [micro sign]g/kg was given, PE^CO(2) exceeded 65 mmHg, or apneic periods lasting more than 60 s occurred. Based on pharmacokinetic simulations, the cumulative dose was chosen to yield concentrations between 200 and 250 [micro sign]g/l, which were identified as allowing adequate spontaneous ventilation during recovery from anesthesia. [15] During the study, an observational scale (1 = awake, restless; 2 = awake, calm; 3 = lightly sedated; 4 = asleep) and continuous electroencephalographic monitoring were used to assess the patients' state of arousal. The study was discontinued at any time when verbal stimulation was required to maintain spontaneous ventilation or the patients were uncomfortable while breathing through the mask. Otherwise data were collected for 60 min after the infusion began. Anesthesia was induced immediately after the study was discontinued.

Arterial blood samples were drawn before and 3, 6, 9, 12, 15, 20, 25, 30, 35, 40, 45, 50, and 60 min after the start of the infusion for alfentanil assay and determination of the Pa^CO(2).

All samples were stored on ice after collection. Before the study, we ensured that storing the arterial blood samples on ice for as long as 2 h in the tubes used for the study did not lead to appreciable changes of the Pa (^CO)(2). We did this in response to several published reports regarding this issue. [16–20] We performed blood gas analyses immediately after the study was discontinued using an automated analyzer with an autocalibration function (ABL 505; Radiometer Medical A/S, Copenhagen, Denmark). The blood samples drawn to analyze the alfentanil concentration were centrifuged at 3,000 rpm for 15 min and the plasma was stored at -20 [degree sign]C until the assay.

Alfentanil concentrations were determined using a sensitive and selective high-performance liquid chromatography assay. Briefly, alfentanil and sufentanil (Janssen, Beerse, Belgium), the latter serving as an internal standard, were extracted from the plasma samples as previously described. [21] The recovery rate was 89%. The samples were analyzed by high-performance liquid chromatography (mobile phase: 0.03 M NaH²PO⁴, 55%; acetonitrile, 35%; isopropanol, 10%; column: Supelcosil DB-8 [250 x 4.6 mm] with ultraviolet detection [220 nm]). Retention times were 5.4 min for alfentanil and 8.8 min for sufentanil. The limit of quantitation was 10 [micro sign]g/l. The assay was linear in a concentration range of 10–1,000 [micro sign]g/l (coefficient of correlation = 0.991). Precision of the assay was 2.95%(coefficient of variation).

The program system NONMEM, version IV with the First Order Conditional Estimation method and [Greek small letter eta]-[Greek small letter epsilon] interaction to reduce the influence of model misspecification was used for all model fits, empirical Bayesian estimation of the individual parameters, and simulations. [22]

The pharmacokinetic-dynamic analysis was performed sequentially. The population means, the interindividual variability, and the empirical Bayesian estimates of the individual pharmacokinetic parameters were determined first. Subsequently, the pharmacodynamic data (Pa^CO(2)) were analyzed with fixed individual pharmacokinetic parameters.

An exponential model was used to describe the interindividual variability in the pharmacokinetic and the pharmacodynamic parameters:Equation 1where [Greek small letter theta] sub (n,i) refers to the individual value of the nth parameter in the ith individual, [Greek small letter theta] sub (n,m) is the typical value in the population of the nth parameter, and [Greek small letter eta] varies randomly among individuals with a mean of zero and a diagonal variance-covariance matrix Omega².

A multiplicative (constant CV) error model was chosen to model residual variability:Equation 2DV^obsrefers to the observed value of the dependent variable (alfentanil concentration, Pa^CO(2)), and DV^exprefers to the value predicted based on dose, time, and the individual pharmacokinetic and pharmacodynamic parameters. [Greek small letter epsilon] is a normally distributed random variable with mean zero and variance [Greek small letter sigma]².

Decisions between different models were made using the Akaike information criterion. [23] We checked for model misspecification by plotting the ratio of the measured and the predicted concentrations against time on a semilogarithmic scale. The covariates tested were weight, age, and completion of the infusion. Covariates were included in the model if the inclusion significantly improved the log likelihood criterion (P < 0.01). [22]

One-, two-, and three-compartment models were fitted and compared as described before. The models were parameterized in terms of the volumes of distribution, the elimination clearance, and the intercompartmental (distribution) clearances.

Because of the direct relation between the volume of a gas and its molar weight (22.4 1 of an ideal gas [22.2 l carbon dioxide] equals 1 mol during standard temperature pressure dry conditions [24]), mass balance equations can be used to describe volume changes. Constant temperature, volumes, and pressures are directly proportional; therefore, changes of partial pressures of a gas also can be computed by mass balance equations. The most simplistic model of carbon dioxide turnover in the body would be a one-compartment model with constant input (carbon dioxide production) and constant output (carbon dioxide elimination) during baseline steady state conditions. The "concentration" in the compartment equals the Pa^CO(2) normalized to atmospheric pressure in our model. Equation 3where V^dCO²is apparent volume of distribution of carbon dioxide in the body (in l); Pa^CO(2) is arterial partial pressure of carbon dioxide (in mmHg); k (ⁱⁿ) is production rate of carbon dioxide (in l/min); k out is elimination rate of carbon dioxide (in l/min); 760 is atmospheric pressure at sea level (in mmHg); k^out(t) is the product of the (pulmonary) clearance (Cl; in l/min) and the Pa^CO(2) at time t divided by the barometric pressure. Equation 4

Considering that the difference between PE^CO(2) and alveolar carbon dioxide pressure (PA^CO(2)) is small in healthy persons, clearance of carbon dioxide can be substituted by alveolar ventilation (V^alv). Equation 5

During baseline steady state conditions (at t = 0), before drug exposure:Equation 6

And therefore Equation 7

Under the assumption that carbon dioxide production remains constant, kⁱⁿ(t) can be expressed as the product of the baseline value of the normalized Pa^CO(2) and alveolar ventilation. Equation 8

Rearranging for the change of Pa^CO(2), the dependent variable, over time yields Equation 9

Alternatively, the ratio of the alveolar ventilation and the apparent volume of distribution of carbon dioxide can be expressed as the elimination constant of carbon dioxide, k^el(min^-1). Equation 10

The Equation mustbe completed for the respiratory depressant effect of opioids and the respiratory stimulant effect of carbon dioxide. Because the volume of distribution remains constant, all changes in carbon dioxide elimination must be caused by changes of V^alv. Opioids are known to reduce (alveolar) ventilation. [25] Because the change in Pa^CO(2) is not a direct effect of opioids but rather is an indirect effect after the reduction of alveolar ventilation caused by opioids, we chose to model the relation between alfentanil concentration and Pa^CO(2) using an indirect-response model. [13,14] A linear model with a negative slope and fractional E^maxmodel were tested. Equation 11, Equation 12with c referring to the current opioid concentration, n referring to the slope factor, and EC⁵⁰referring to the concentration at which V^alv, and therefore k^el, would be decreased to 50% of the value in the absence of opioid. To avoid redundancy, the further derivations are shown for the fractional E^maxmodel only.

Carbon dioxide is a powerful respiratory stimulant. [1] Its response curves are shaped like a hockey stick, [1,24] which is usually described in terms of a carbon dioxide threshold and the linear slope of the carbon dioxide response curve. The following hyperbolic function displays a similar shape and predicts a continuous, nonlinear decrease of the ventilation at less than physiologic Pa^CO(2) which corresponds with a study that evaluated the respiratory rate and inspiratory effort during hypocapnic and hyperoxic conditions in awake volunteers. [26]Equation 13

F denotes the gain of the system. Equation 12 and Equation 13can be combined to describe the net effect for any given opioid concentration and Pa^CO(2) on alveolar ventilation as the product of the terms. This is depicted in Figure 1. Equation 14

Note that an opioid concentration equal to the EC⁵⁰only leads to a 50% reduction in V^alvif the output of the second term equals 1.

Combining Equation 14and Equation 19 yields the final Equation todescribe opioid-induced hypercapnia. Equation 15

A combined pharmacokinetic-pharmacodynamic simulation of the hypercapnic effect of alfentanil was performed for a population of 100 participants during an intravenous infusion of 5.78 mg in 30 min, the calculated dose for the median participant (Wt = 82.5 kg). One hundred sets of individual pharmacokinetic and pharmacodynamic parameters were simulated based on the estimated population means and interindividual variances.

Further simulations were performed to obtain estimates of the steady state Pa^CO(2), corresponding to several alfentanil concentrations (50, 100, 150, 200, and 250 [micro sign]g/l). Baseline alveolar ventilation was not measured, so we cannot calculate absolute alveolar ventilation. However, because any change of k^elmust be attributed to a change in V (^alv), the relative change of alveolar ventilation can be calculated as shown below to obtain steady state alveolar ventilation normalized to baseline during different opioid concentrations. Equation 16The subscript "ss" refers to steady state. Alveolar ventilation normalized to baseline is now called "fractional alveolar ventilation."

In 6 of 14 patients studied, the infusion had to be discontinued prematurely. Three experienced apneic periods lasting more than 60 s, and the other three reached an PE^CO(2) more than 65 mmHg. In three patients, the data collection had to be discontinued prematurely: one needed verbal stimulation to continue breathing after 20 min of drug administration and two patients were uncomfortable and nauseated while they breathed through the mask after 45 min. They received 1.25 mg droperidol as an intravenous bolus to treat nausea, both after 46 min. One patient received 10 mg intravenous urapidil, a substance with [Greek small letter alpha]^{1-antagonistic}and [Greek small letter alpha]^2-agonisticproperties, to treat hypertension after 17 min. No patient was more than lightly sedated, as judged from clinical observation and electroencephalographic data.

The concentration time course of alfentanil in 60 min was adequately described by a two-compartment model. Table 1shows the pharmacokinetic parameters. Figure 2, A and B show the goodness of fit under the population model and using the Bayesian estimates of the individual pharmacokinetic parameters. As shown in Figure 2C, the population means and variability parameters of the mixed-effects model were used to simulate a population of 100 persons receiving the median cumulative dose (5.78 mg in 30 min). The actual measured concentration-time data are superimposed on the simulations. The data from persons who did not receive the full dose are included up to the discontinuation of the infusion. We must stress that all concentration-time courses are predicted well with Bayesian estimates of the pharmacokinetic parameters (Figure 2B), because these were used to generate the concentration-time profiles that provided the input for the pharmacodynamic model.

(Figure 3) depicts the relation between plasma concentrations predicted by dose, time, and the Bayesian estimates of the individual pharmacokinetic parameters and the measured Pa^CO(2) concentrations. A counterclockwise hysteresis can be observed.

Both by inspection of the residuals and by the Akaike information criterion, the fractional E^maxmodel was found to be superior to the negative-slope model for describing the respiratory depressant effect of alfentanil. Table 2shows the pharmacodynamic parameters of the E^maxmodel. Because we could not simultaneously estimate F and EC⁵⁰, the population mean of F was fixed, as indicated in Materials and Methods. The interindividual variability of F was estimated using the model. There was a trend toward lower values of F in the post hoc estimates of the six patients in whom the infusion was discontinued prematurely, which was not statistically significant. Figure 4A shows the goodness of fit under the population model, and Figure 4B shows the goodness of fit using the Bayesian estimates of the individual pharmacodynamic parameters. In Figure 4C, the population means and variability parameters of the mixed-effects model were used to simulate a population of 100 persons who received the median cumulative dose (5.78 mg in 30 min). The measured Pa^CO(2) data are superimposed on the simulation. The data from persons who did not receive the full dose are included until the infusion was stopped.

(Figure 5A) shows the predicted steady state relation between alfentanil concentration and Pa^CO(2), and Figure 5B shows the corresponding fractional alveolar ventilation calculated as indicated in Materials and Methods. One hundred subjects per concentration were simulated. Displayed are the individual values, the median, and the 16 and 84% quartiles of the distribution for steady state concentrations within the range of the concentrations observed during the study. The median values are also summarized in Table 3.

(Figure 6, A and B) show the relation of fractional alveolar ventilation, Pa^CO(2), and opioid concentration for the mean participant. In Figure 6A, the mean response surface, all data points obtained (dots), and the calculated mean steady state ventilation (bold line) are displayed. In Figure 6B, the mean response surface has been shaded and the individual data and steady state response are omitted.

The goal of this study was to provide a model-based description of the time course of opioid-induced respiratory depression. Therefore, we did not obtain blood samples after the pharmacodynamic measurements were made. The pharmacokinetic model merely serves as the input for the pharmacodynamic model during the Pa^CO(2) measurements. The parameters cannot be used to extrapolate the concentration time course of alfentanil beyond the study period.

The following aspects must be considered when modeling opioid-induced respiratory depression: the end point assessed, the design of the study, and the modeling approach used. The effects of opioids on ventilation are classically described using carbon dioxide response curves. [1–3] We decided against this approach for three reasons. First, because recording a single carbon dioxide response curve takes minutes, the drug concentration at the effect site must be held steady. This can be achieved with computer-controlled drug delivery but limits the number of concentrations that can be assessed. Second, carbon dioxide response curves require the administration or rebreathing of carbon dioxide, which must be allowed to wash out before a new measurement can be obtained. By allowing one measurement every 15 min, [3] this severely limited the resolution of the method. Although this design was used by Hill et al. [27] to compare the respiratory depressant action of opioids with considerable success, it does not provide sufficient points on the concentration-effect curve for a modeling approach. Furthermore, steady state data are not helpful when trying to characterize the behavior of the system in a non-steady state situation, which is commonly encountered in the clinical setting. For the latter reason, our goal was to perform a non-steady state analysis. Third, the clinical significance of the parameters obtained is difficult to assess. Clinicians will be considerably more interested in the opioid concentration that increases the Pa^CO(2) or that decreases alveolar ventilation to a certain value than in the concentration that decreases the slope of the carbon dioxide response curve by some percentage. Exactly the same problems will be encountered when measuring the opioid-induced attenuation of the hypoxic respiratory drive. For these reasons, carbon dioxide-oxygen response curves were excluded a priori as a primary end point. Minute ventilation is influenced by Pa^CO(2) and opioid concentration, as can be seen from carbon dioxide response curves under different opioid concentrations. The same is true for alveolar ventilation. Therefore, modeling the influence of an opioid on minute or alveolar ventilation without concomitantly measuring and modeling Pa^CO(2) or PE^CO(2) will lead to a biased estimate for the potency of the opioid. Because the respiratory stimulant effect of carbon dioxide would be missed, falsely high estimates of the EC⁵⁰would be obtained. Even worse, because different input functions of the opioid, most notably bolus versus slow infusion, would lead to different time courses of Pa^CO(2), which would not be accounted for in the model, the estimate of the EC⁵⁰might become dependent on the experimental design. These systematic errors might occur even in the presence of an adequate fit for the data measured.

The Pa^CO(2) is easily measured; the measurements can be taken in short intervals (minutes); and it is highly clinically significant. Therefore, we chose Pa^CO(2) as our pharmacodynamic end point. To achieve a gradually increasing degree of respiratory depression and pronounced hypercapnia without encountering respiratory arrest, we decided to administer the drug by zero-order infusion. We aimed for a target concentration between 200 and 250 [micro sign]g/l, which was previously determined to allow adequate spontaneous ventilation after surgery. [15] Although the infusion was discontinued prematurely in 6 of 14 patients, only one patient was prematurely lost to data acquisition because of pronounced respiratory depression. Conversely, nearly every patient showed pronounced hypercapnia (Figure 4C). We believe the goal of obtaining a sufficient range of the Pa (^CO)(2) without compromising patient safety has been adequately achieved.

As stated in Materials and Methods, the relation between alfentanil concentration and Pa^CO(2) is indirect. As expected, a modified indirect-response model adequately described the data. Although effect compartment concentrations had to be calculated to model the electroencephalographic slowing effect of alfentanil, [6] alfentanil-induced respiratory depression could be modeled assuming a direct relation between plasma concentrations and effect. Because the plasma effect compartment equilibration half-time of alfentanil is 1.1 min, much smaller than the 5 min required for a 50% change in the Pa^CO(2) from the current to the new steady state level after a change of alveolar ventilation, [28] the k^eoof alfentanil could not be estimated from our data. The hysteresis displayed in Figure 3was explained entirely by the inherent inertia of the indirect effect model. However, failure to account for the equilibration of alfentanil with the effect site could have distorted the estimate of the elimination constant of carbon dioxide (k^el). Because the product of k^eland the baseline value of Pa^CO(2)(Pa^CO(2)(0)), yields the fastest possible increase of Pa^CO(2) according to the model, it can be used to check for plausibility. The model predicts a maximal rate of 3.6 mmHg/min, which agrees closely with the value published in a standard text of respiratory physiology (3–6 mmHg/min). [24]

Because we could assign a physiologic meaning to the elimination constant of the indirect response model (the ratio of alveolar ventilation and the volume of distribution of carbon dioxide), we could predict changes in alveolar ventilation using parameters estimated solely from fitting Pa^CO(2) data. Because patients with opioid-induced respiratory depression are more likely to have complications from hypoxemia than from hypercapnia, this feature of the model is especially attractive.

The following potential shortcomings of our approach must be considered. First, as can be seen in the derivation of the model, the production rate of carbon dioxide has been substituted by the product of the elimination rate and Pa^CO(2) in the absence of drug, which is equal to the assumption that carbon dioxide production remains constant during the entire study. In a recent study, carbon dioxide production was diminished 12–19% from baseline during a target-controlled infusion of alfentanil to 40 [micro sign]g/l, [29] a very low concentration compared with the 300 [micro sign]g/l observed in some patients in our study. However, there is clear evidence that this trend does not continue. Hill et al. [27] measured carbon dioxide production of healthy volunteers at a steady state alfentanil concentration of 80 [micro sign]g/l and found no significant difference from baseline. Weyland et al. [30] compared oxygen consumption during sedation with midazolam, propofol, thiopental, and fentanyl. All hypnotics reduced oxygen consumption by as much as 15% from baseline. Patients receiving fentanyl exhibited a slight (5%) but not significant increase in oxygen consumption, at a Pe^CO(2) of 50 mmHg, when compared with placebo. Interestingly, these patients were only lightly sedated at this degree of respiratory depression, which corroborates our findings that profound opioid-induced respiratory depression is not accompanied by substantial sedation or sleep.

Assuming a constant carbon dioxide production during the course of the study cannot be considered a relevant systematic error. Nevertheless, by introducing several pseudo-steady state steps with a computer-controlled infusion and measuring carbon dioxide production before drug administration and during pseudo-steady state, this potential pitfall could be avoided completely.

Second, because we could not simultaneously fit the population means of both the EC⁵⁰and F, the gain of the carbon dioxide response curve, we had to fix the mean response to carbon dioxide in the population to a predetermined value. Because we estimated the interindividual variability in the population, we were able to individualize the gain of the carbon dioxide response curve. This technique, assigning a mean response and estimating the magnitude of the interindividual variability, can be used only with a population approach. However, measurement of a single carbon dioxide response curve in the absence of opioid in every patient before the study would have directly provided information about the individual sensitivity to carbon dioxide, which could have been exploited in the modeling approach. This fact should be accounted for in future designs.

Third, our modeling approach considers the body as a single compartment, with a partial pressure of carbon dioxide equaling Pa^CO(2). Because it is obvious that the Pa^CO(2) and central venous partial pressures of carbon dioxide differ, this simplification contradicts existing knowledge of physiology. The derivation of a two-compartment model to acknowledge this fact can be achieved easily. Because we do not have any evidence of model misspecification (Figure 4), the simple model accounts for the observed time course of carbon dioxide, which means that the information content of the data is not sufficient to estimate a higher number of parameters than used in a one-compartment model. Therefore, we favor the more parsimonious approach.

With regard to the interpretation of the results obtained, the following points should be made. First, as can be seen in Figure 5, relatively high plasma concentrations of alfentanil can be tolerated during steady state conditions. The corresponding fractional alveolar ventilation of 0.72 (median) for a concentration of 250 [micro sign]g/l clearly provides adequate oxygenation with an FI^O(2) of 0.21 and atmospheric pressure at sea level (alveolar oxygen tension = approximately 82 mmHg, according to the alveolar gas equation [24]). However, the fractional alveolar ventilation of several persons at that concentration is closer to 0.5 (alveolar oxygen tension = approximately 66 mmHg), a clearly undesirable degree of respiratory depression. Furthermore, considering that similar concentrations were to be reached instantaneously, a precipitous decrease in the alveolar ventilation and, as a logical consequence, hypoxia, would certainly occur, as can be deduced from Figure 6. Immediately after an alfentanil (opioid) bolus, high concentrations of the drug at a normal Pa^CO(2) would lead to a pronounced decrease of alveolar ventilation. Thereafter, because carbon dioxide accumulation occurs and is a result of its respiratory stimulant effect, alveolar ventilation would increase again, even if the alfentanil (opioid) concentration had been kept constant. The arrows in the Figure areintended to portray this situation.

Can we predict a safe alfentanil concentration for the spontaneously breathing patient? With the caveats of not taking into account anesthetic drug interactions and the stimulatory effect of pain on respiration, [31–33] a steady state concentration of 50 [micro sign]g/l will decrease alveolar ventilation less than 23% from baseline in 97.5% of the patients. Steady state Pa^CO(2) would be less than 53.6 mmHg in 97.5% of the patients. Regarding the fact that reported minimal effective concentrations for postoperative analgesia range from 15 [micro sign]g/l to 60 [micro sign]g/l, [34,35] this recommendation looks equally attractive when based on effectiveness data obtained in patients after surgery.

In conclusion, we developed a pharmacokinetic-pharmacodynamic model for the quantitative description of opioid-induced respiratory depression. As noted previously, several important questions must be considered regarding opioid-induced respiratory depression, and they could be answered quantitatively by applying the methods we described.

The authors thank Dr. Elizabeth Youngs and Dr. Steven Shafer, Stanford University, for critical discussions of the manuscript.