To simulate the time course of drug effect, it is sometimes necessary to combine the pharmacodynamic parameters from an integrated pharmacodynamic-pharmacodynamic study (e.g., volumes, clearances, k(e0) [the effect site equilibration rate constant], C(50) [the steady state plasma concentration associated with 50% maximum effect], and the Hill coefficient) with pharmacokinetic parameters from a different study (e.g., a study examining a different age group or sampling over longer periods of time). Pharmacokinetic-pharmacodynamic parameters form an interlocked vector that describes the relationship between input (dose) and output (effect). Unintended consequences may result if individual elements of this vector (e.g., k(e0)) are combined with pharmacokinetic parameters from a different study. The authors propose an alternative methodology to rationally combine the results of separate pharmacokinetic and pharmacodynamic studies, based on t(peak), the time of peak effect after bolus injection.

The naive approach to combining separate pharmacokinetic and pharmacodynamic studies is to simply take the k(e0) from the pharmacodynamic study and apply it naively to the pharmacokinetic study of interest. In the t(peak) approach, k(e0) is recalculated using the pharmacokinetics of interest to yield the correct time of peak effect. The authors proposed that the t(peak) method would yield better predictions of the time course of drug effect than the naive approach. They tested this hypothesis in three simulations: thiopental, remifentanil, and propofol.

In each set of simulations, the t(peak) method better approximated the postulated "true" time course of drug effect than the naive method.

T(peak) is a useful pharmacodynamic parameter and can be used to link separate pharmacokinetic and pharmacodynamic studies. This addresses a common difficulty in clinical pharmacology simulation and control problems, where there is usually a wide choice of pharmacokinetic models but only one or two published pharmacokinetic-pharmacodynamic models. The results will be immediately applicable to target-controlled anesthetic infusion systems, where linkage of separate pharmacokinetic and pharmacodynamic parameters into a single model is inherent in several target-controlled infusion designs.

INTEGRATED pharmacokinetic–pharmacodynamic studies require carefully gathering blood samples, making frequent measurements of drug effect, and developing complex models of the full pharmacokinetic and pharmacodynamic behavior of the drug. 1Because integrated pharmacokinetic–pharmacodynamic studies are so difficult and because drug assays typically improve considerably over the useful life of a new drug, the typical research history for an intravenous anesthetic drug involves many pharmacokinetic studies but only one or two integrated pharmacokinetic–pharmacodynamic studies. As a result, the definitive pharmacokinetic study in a population of interest may not incorporate a pharmacodynamic component (*e.g.* , sufentanil pharmacokinetics reported by Gepts *et al.* , 2propofol pharmacokinetics in children reported by Kataria *et al.* , 3propofol pharmacokinetics during cardiopulmonary bypass reported by Bailey *et al.* 4).

Effective use of “pharmacokinetic-only” studies in simulations requires linking them to pharmacodynamic models. The naive way to perform this linking is to simply take the value of k^{e0}, the plasma effect site equilibration rate constant, reported in the integrated pharmacokinetic–pharmacodynamic study, and use this value with the pharmacokinetic parameters of interest to simulating the time course of drug concentration in the effect site. As Gentry *et al.* 5demonstrated for thiopental, and Wakeling *et al.* 6demonstrated for propofol, this naive approach is unwise and will lead to very poor predictions of the time course of drug effect.

We propose an alternative approach for combining pharmacodynamic parameters from an integrated pharmacokinetic–pharmacodynamic study with pharmacokinetic parameters from a separate study, based on the time of peak effect, t^{peak}. In this approach, the investigator first translates the k^{e0}reported in the integrated pharmacokinetic–pharmacodynamic study into t^{peak}, the time of maximum effect site concentration following an intravenous bolus dose when there is no drug initially in the system. Of course, if the original study observes and reports t^{peak}, this can be used instead. The investigator then calculates the value of k^{e0}that accurately predicts t^{peak}when used with the pharmacokinetic parameter set of interest. We propose that t^{peak}is a model independent parameter, as it can be directly observed following a submaximal bolus dose. As a model independent parameter, t^{peak}can be used with different pharmacokinetic parameter sets to accurately describe the time course of concentration in the effect site.

The goal of this article is to evaluate whether the proposed “t^{peak}method” is superior to the naive approach with three often used drugs, remifentanil, thiopental, and propofol, as follows:

To compare the predicted effect site concentrations obtained by the two methods, we evaluate two previously published studies of thiopental. 7,8Each of these studies is an integrated pharmacokinetic–pharmacodynamic study, which allowed examination of the results of exchanging parameters between the studies.

We performed a simulation study based on the pharmacokinetics and pharmacodynamics of remifentanil 9to compare the time course of drug effect and the accuracy of k

^{e0}estimates for the two methods.We applied the t

^{peak}method to propofol to combine propofol pharmacokinetic parameters 10that have been demonstrated to be reliable 11with the k^{e0}from a published pharmacokinetic–pharmacodynamic model. 12,13These results are directly relevant to the commercially available propofol infusion system Diprifusor (

AstraZeneca, London, United Kingdom).

## Methods

As a theoretical treatise, this work was not reviewed by an institutional review board. The published studies explored in this analysis were approved by their respective institutional review boards. 7–13

We suppose that we have available the results from the “best” combined pharmacokinetic–pharmacodynamic study and have selected a “best” pharmacokinetic study that addresses our specific population. The proposed method for combining the results of these two studies is:

From the parameters of the “best” combined pharmacokinetic–pharmacodynamic study, calculate t

^{peak}, the time of peak effect following submaximal bolus administration. Of course, if t^{peak}has been observed and reported as part of the original study, that value can be used directly.Using the “best” pharmacokinetic parameters calculate an “updated” k

^{e0}that preserves the t^{peak}calculated in step 1 or reported in the original study.

The mathematical details of calculating t^{peak}are explored in the Appendix. An intuitive representation follows. Let the concentration of drug in the plasma following a bolus dose of 1 unit have the polyexponential form

with the bolus occurring at time 0, when there is no drug in the system. Note that this model assumes instantaneous mixing. The limitations of this assumption are presented in the Discussion.

Let the concentration of drug in the effect site be related to the concentration in the plasma by a first-order process with rate constant k^{e0}. Note that this assumes a “direct effect” model. Again, the limitations of this assumption are presented in the discussion.

With these assumptions, the concentration in the effect site resulting from a bolus dose of 1 unit is given by The time t^{peak}is the time at which C^{e}reaches its maximum value. Since we are only concerned with the time of maximum effect site concentration, the magnitude of the bolus dose is irrelevant: changing the magnitude changes the maximum concentration, but not the time at which it occurs. Thus, we can think of C^{p}(t) as the unit disposition function (UDF) of the plasma, and we can think of C^{e}(t) as the UDF of the effect site. To use t^{peak}in the manner we have described, we are presented with two problems: finding t^{peak}given k^{e0}, and finding k^{e0}given t^{peak}.

In the first problem, finding t^{peak}means finding the time t that maximizes C^{e}in equation 2, given the pharmacokinetic parameters A and λ. This maximization problem can be solved numerically. For the second problem, finding k^{e0}given t^{peak}and a new vector of pharmacokinetic parameters A and λ, we observe that since t^{peak}maximizes C^{e}, the derivative of C^{e}with respect to t must equal zero at t = t^{peak}. That is, differentiating equation 2, t^{peak}satisfies

Thus, finding k^{e0}means solving equation 3for k^{e0}, which is done numerically. The Appendixcovers a number of important technical mathematical details of finding t^{peak}given k^{e0}and a vector of pharmacokinetic parameters, and finding k^{e0}given t^{peak}and a different vector of pharmacokinetic parameters.

### Simulation studies

#### Thiopental.

To compare the naive method with the t^{peak}method of combining pharmacokinetic data from one study with pharmacodynamic data from a second study, we considered the models developed in two previously published studies of thiopental, one by Stanski and Maitre 7and the other by Shanks *et al.* 8The two studies are summarized in table 1. In the model of Stanski and Maitre, 7we chose parameters corresponding to a 35-yr-old, 70-kg subject. The UDF (thiopental concentration, in units of μg/ml, in the plasma following a bolus dose of 1 mg) corresponding to these values is

In the model by Shanks, 8the plasma thiopental concentrations over time are described by the following UDF:

Based on these pharmacokinetic models, we constructed two hypothetical cases and evaluated each using the naive and t^{peak}methodology:

### Test Case 1

We imagine that we have the pharmacokinetic and pharmacodynamic results from Shanks *et al.* , 8and we wish to use their pharmacodynamics (*i.e.* , k^{e0}or t^{peak}) with the “improved” age- and weight-adjusted pharmacokinetic model of Stanski and Maitre. 7In this test, the full pharmacokinetic–pharmacodynamic model reported by Stanski and Maitre 7represents the “hypothetical truth,” and the test is whether the naive approach or the t^{peak}approach better integrates the pharmacodynamics of Shanks *et al.* 8with the pharmacokinetics of Stanski and Maitre. 7

### Test Case 2

We imagine that we have the pharmacokinetic and pharmacodynamic results from Stanski and Maitre, 7and we wish to use their k^{e0}with the “improved” four-compartment pharmacokinetic model of Shanks *et al.* 8In this test, the full pharmacokinetic–pharmacodynamic model reported by Shanks *et al.* 8represents the “hypothetical truth,” and the test is whether the naive approach or the t^{peak}approach better integrates the pharmacodynamics of Stanski and Maitre 7with the pharmacokinetics Shanks *et al.* 8

For both test cases, the goodness of the naive and t^{peak}approaches was assessed graphically.

### Remifentanil

We used simulations of remifentanil, for which we have detailed population models, 9to address the accuracy of the naive and t^{peak}methods of determining the k^{e0}when an improved pharmacokinetic model is developed. We created 100 individuals drawn from the population distribution as described by Minto *et al.* 9and defined as “true” the time course of plasma and effect site concentration following a 500-μg intravenous bolus of remifentanil. We then simulated a crude pharmacokinetic–pharmacodynamic clinical trial as follows:

Simulated blood samples were drawn at 0, 5, 10, 15, 20, 25, 30, 40, and 60 min. No noise was added to the assay.

Drug effect was measured concurrent with the blood samples. No noise was added to the “measured” drug effect.

A one-compartment pharmacokinetic model was fit to the measured remifentanil concentrations.

A k

^{e0}link model and a sigmoid E^{max}pharmacodynamic model were fit to the measured drug effect.Two determinations of t

^{peak}were made:a. Calculated from the model built in step

b. Directly “observed” during the experiment

We then assumed that an improved pharmacokinetic model was developed in a subsequent pharmacokinetic study. In this case, the improved model was simply the “true” remifentanil pharmacokinetics for each individual, as described above. We examined the time course of drug effect when the true pharmacokinetic model was combined with three different estimates of k^{e0}:

Naive method: k

^{e0}determined from the crude studyt

^{peak}method: t^{peak}calculated from the crude study (4a above)t

^{peak}method: t^{peak}directly observed (4b above)

Goodness of each method was assessed graphically and by repeated-measures *t* test of the k^{e0}value *versus* the true value of k^{e0}in each subject.

### Propofol

The third test has immediate commercial implications: how to best combine the propofol pharmacokinetic parameters reported by Marsh *et al.* , 10as presently used in the Diprifusor device, with the pharmacodynamic parameters from the pharmacokinetic–pharmacodynamic model published by Schnider *et al.* 12,13The pharmacokinetic model reported by Marsh *et al.* 10has been demonstrated to produce acceptable performance when used to control the administration of propofol 11and thus is the pharmacokinetic model incorporated into the only commercially approved (outside of the United States) propofol target-controlled infusion pump (Diprifusor). Marsh *et al.* 10did not concurrently estimate pharmacodynamic parameters, and thus, the question of how to combine the Marsh pharmacokinetic parameters with a propofol pharmacodynamic model is of immediate consequence to the Diprifusor device.

From the pharmacokinetic model reported by Schnider *et al.* , 12we chose parameters corresponding to a 35-yr-old, 70-kg, 178-cm man. The UDF (propofol concentration, in units of μg/ml, in the plasma following a bolus dose of 1 mg) corresponding to these values is

From the model reported by Marsh *et al.* , 10we chose parameters corresponding to a 70-kg subject. The UDF corresponding to these values is

The question is: What is the best value of k^{e0}to use with the Marsh pharmacokinetics? We assumed, somewhat parochially, that the high-resolution pharmacokinetic–pharmacodynamic study of Schnider correctly characterized the time course of effect site propofol concentration following bolus injection. We compared two different methods of deriving k^{e0}for use with the Marsh pharmacokinetics:

Naive method: The k

^{e0}reported by Schnider 13is used with the pharmacokinetics reported by Marsh. 10t

^{peak}method: k^{e0}is calculated for the Marsh pharmacokinetics 10so that the peak effect site concentration occurs at the same time as the peak effect site concentration reported by Schnider. 13

The goodness of each method was assessed graphically.

All modeling was done with Excel (Microsoft, Redmond, WA), using the “Solver” function and macro routines written by the authors.

## Results

### Thiopental Simulation Study

Table 1summarizes the thiopental studies used in the simulations. Figure 1shows the results of test case 1, in which the “truth” was defined as the pharmacokinetic and pharmacodynamic results from Shanks 8(solid line). The dotted line shows the results of the naive approach, in which the k^{e0}from the Stanski and Maitre 7trial (0.58 min^{−1}) has been combined with the pharmacokinetics reported by Shanks *et al.* 8The naive approach poorly approximates the “true” time course of effect site concentrations. The dashed line shows the result of the t^{peak}approach. The Stanski and Maitre 7pharmacokinetic–pharmacodynamic model predicts a t^{peak}of 1.75 min. This t^{peak}, when combined with the pharmacokinetic parameters of Shanks *et al.* , 8yielded a k^{e0}of 0.284 min^{−1}, in close agreement with the value of 0.29 min^{−1}reported by Shanks *et al.* 8The dashed line is nearly indistinguishable from the solid line that defines the true model, demonstrating that the t^{peak}approach has almost perfectly captured the “true” time course of drug effect.

Figure 2shows the results of test case 2, in which the “truth” was defined at the pharmacokinetic and pharmacodynamic results as reported by Stanski and Maitre 7(solid line). The dotted line shows the results of the naive approach, in which the k^{e0}reported by Shanks *et al.* 8(0.29 min^{−1}) has been combined with the pharmacokinetics reported by Stanski and Maitre. 7The naive approach poorly approximates the “true” time course of effect site concentration. As in figure 1, a dashed line in figure 2shows the prediction using the t^{peak}approach, in which the time of peak effect reported by Shanks *et al.* 8has been used to calculate a k^{e0}for use with the pharmacokinetics reported by Stanski and Maitre. 7The dashed line is not seen in figure 2because it exactly matches the “true” model, and thus is hidden under the solid line.

In this example, the t^{peak}approach accurately reproduced the “true” time course of drug effect in both simulations. The naive approach did not accurately reproduce the true time course in either simulation.

### Remifentanil Simulation

The top graph in figure 3shows the true values of the plasma concentration (dotted lines) and the effect site concentration (solid lines) for the 100 simulated individuals. The drug effect peaks approximately 1–2 min after bolus injection. The level of interindividual variability is small, reflecting the low variability in remifentanil pharmacokinetics. 9The average k^{e0}was 0.69 ± 0.039 (SEM) min^{−1}. The lower graph in figure 3shows the result of the poor pharmacokinetic–pharmacodynamic study in the same 100 individuals. Because the first sample after baseline was taken at 5 min, the initial high remifentanil concentrations were not captured in the pharmacokinetic model, affecting both the apparent time course of remifentanil concentration and the magnitude of the effect site concentrations. The average k^{e0}in the crude study was 2.1 ± 0.15 (SEM) min^{−1}, which was significantly larger than the true values of k^{e0}(*P* << 0.001).

Figure 4shows the results of the three different approaches to determining the right k^{e0}, given the correct pharmacokinetics in each individual but only the pharmacodynamic results from the crude pharmacokinetic–pharmacodynamic trial shown in the lower graph of figure 3. In the naive approach (top graph), the effect site concentrations uniformly equilibrate too quickly with the plasma, resulting in a predicted effect site concentration that rises to a higher, earlier peak than predicted by the “true” model (fig. 3, top graph). This is the expected result of the consistently larger value of k^{e0}in estimated in the crude pharmacokinetic–pharmacodynamic study. The middle graph in figure 4shows the predicted effect site concentrations over time using the t^{peak}approach, in which k^{e0}is calculated for each individual to provide the same time of peak effect as predicted by the model parameters of the crude pharmacokinetic–pharmacodynamic study. Despite the crude nature of the study, the time course of concentration in the effect site closely matches the “true” time course shown in the top graph of figure 3. k^{e0}estimated using the calculated t^{peak}was 0.55 ± 0.06 min^{−1}, which was not statistically different from the true values of k^{e0}. The bottom graph in figure 4shows the predicted effect site concentrations over time using the value of k^{e0}that matched the “observed” time of peak effect (*i.e.* , the time of peak effect seen in the top graph of fig. 3). The value of k^{e0}in every individual matched the “true” value of k^{e0}to five significant digits, which was the resolution of the k^{e0}search algorithm. As a result, the time course of drug effect in the lower graph of figure 4exactly matches the time course shown in the top graph of figure 3.

### Propofol Simulation

Figure 5shows the time course of propofol concentration in the effect site as predicted by the Schnider 12,13integrated pharmacokinetic–pharmacodynamic model (solid line). The dotted line shows the results of the naive approach, in which the pharmacokinetics reported by Marsh 10have been combined with the k^{e0}reported by Schnider. 13The effect site concentrations are too low, and the time of peak effect is almost twice the prediction of the Schnider model. 12,13The dashed line shows the results of the t^{peak}approach, in which a k^{e0}was calculated for the pharmacokinetics reported by Marsh 10to give the same time of peak effect as reported by Schnider. 12,13Although this model does not exactly match the “true” (for these purposes) model reported by Schnider, 12,13it comes far closer than the naive approach.

## Discussion

Our simulation analyses extends those of Gentry *et al.* , 5whose work contributed greatly to our thinking. We have shown with three very different simulations that the t^{peak}method of linking disparate pharmacokinetic and pharmacodynamic models is superior to the alternative of simply combining the k^{e0}from one study with the pharmacokinetic parameters from another study. In every simulation, the t^{peak}method yielded k^{e0}estimates closer to the “true” (simulated) k^{e0}values than those obtained by the naive method. As expected, this resulted in a time course of drug effect that also matched the “true” time course of drug effect in each simulation. In the thiopental simulations, the match was so good that the line for the t^{peak}method almost disappeared into the “true” line in the first simulation (fig. 1) and completely disappeared under the “true” line in the second simulation (fig. 2).

None of this should be surprising. It is almost a tautology that linking parameters from two models necessarily will perform worse than an approach that adjusts the parameters to account for differences in the underlying models. This is particularly the case for the time course of drug effect because the estimate of k^{e0}is highly dependent on an accurate description of the pharmacokinetic behavior of the drug. The pharmacokinetic model in turn requires a sound study design and an accurate pharmacokinetic model. If the pharmacokinetic model is “wrong” because of a poor study design or a poor choice of model (*e.g.* , missing the early distribution phase), the estimate of k^{e0}will be affected. In contrast t^{peak}is potentially directly observable following a submaximal intravenous bolus dose. As shown in the third remifentanil example, a directly observed t^{peak}is a model-independent gold standard for describing with a single parameter the time course of concentration in the effect site.

In the remifentanil example, it may seem tautologic that last approach, using the directly observed time of peak effect, would do the best since the time of peak effect was directly observed and not muddled by the poor determination of remifentanil pharmacokinetics. That is exactly the point: As a directly observed pharmacodynamic variable, the time of peak effect should be accurate and independent of the underlying pharmacokinetic and pharmacodynamic model. The remifentanil simulation was structured to make this point obvious.

The choice of the t^{peak}method is not dependent on the structure of the pharmacokinetic model. Mammillary models have well identified problems in the first few minutes, which are addressed by the use of recirculatory models. 14The recirculatory model has a very different structure from the simply mammillary model used in our simulations but fundamentally still relates dose to plasma concentration using a unit disposition function. Therefore, time course of concentration in the effect site concentration can still be modeled as the convolution of the plasma concentrations over time against k^{e0}e^{−keot}, the disposition function of the effect site. The t^{peak}method can be used to find the “correct” value of k^{e0}for use with recirculatory models. One simply finds the value of k^{e0}that predicts the correct time of peak effect following bolus injection.

There are two caveats that need to be mentioned. First, our analysis and the derivation in the Appendixare based on the assumption that the pharmacokinetics are described by an effect site disposition curve with a single peak. It is possible for recirculatory models to have multiple peaks, typically separated by 15–30 s. Given that the typical blood–brain equilibration delay is on the order of at least 30–60 s, the concentrations at the effect site are still likely to be unimodal even with a multimodel plasma drug model. Should injection of a bolus be followed by multiple peaks in the effect site, then it would be necessary to designate the time of peak effect as the time of the first, second, or nth, peak effect. To the best of our knowledge, all anesthetic drugs only produce a single peak drug effect following bolus injection, and so the question of how to handle multiple peaks is not clinically relevant for anesthetic drugs.

Second, we have assumed direct effect models, where the drug effect is an instantaneous reflection of the drug concentration in the effect site. For some effects, such as ventilatory depression, it is much more appropriate to use indirect effect models. 15Indirect models usually contain imbedded direct effect models, and the time to peak effect for the imbedded direct effect model is easily calculated. Although the direct drug effect is almost never observed with indirect effect models, the calculated t^{peak}could still be used to integrate a new pharmacokinetic study with an existing indirect effect model. With indirect effect models, we still encourage investigators to record the time of the observed peak effect (*e.g.* , time of peak ventilatory depression), which provides a model-independent measure of the time course of drug effect that all derivative pharmacokinetic–pharmacodynamic models should accurately predict.

The question remains as to whether it is preferable to choose a model that is the hybrid of a pharmacokinetic–pharmacodynamic model and the “best” pharmacokinetic model, or simply use the “best” of the available integrated pharmacokinetic–pharmacodynamic models. The issue is that the “improved” pharmacokinetic model often represents a major improvement, as explored in the remifentanil example, which motivates the search for how best to combine the results from different studies. There may be practical reasons as well. For example, AstraZeneca carefully explored which of three pharmacokinetic models produced the smallest performance errors for the Diprifusor. 11This evaluation was very expensive, and the pharmacokinetic model selected for the Diprifusor has been approved by a variety of regulatory authorities. This regulatory approval means that the model is not likely to change. However, it would be useful for the next generation of Diprifusor to incorporate the concentrations at the site of drug effect. The only option is to combine pharmacokinetics incorporated in the Diprifusor with the best available model of propofol pharmacodynamics.

Struys *et al.* 16have directly compared the naive and t^{peak}approach to mixing propofol pharmacokinetic and pharmacodynamic models for propofol. They administered propofol to 120 patients, targeting the plasma concentration, the effect site using the naive method for determining k^{e0}, or the effect site using the t^{peak}method for determining k^{e0}. They reported that the “biophase model combining the Marsh kinetics and a time to peak effect of 1.6 min accurately predicted the time course of propofol drug effect.” This is a prospective test of the two simulations shown in figure 5, validating in patients the conclusion of our simulation that the t^{peak}approach should work better.

Many individuals have used the programs STANPUMP #and RUGLOOP **for pharmacokinetic and pharmacodynamic simulations. These programs have used the t^{peak}approach for nearly 10 yr.

In summary, we encourage clinical investigators to design studies in which t^{peak}can be directly observed following bolus injection. This requires administration of submaximal doses so that the peak is not obscured by the drug effect maintaining a plateau at maximum effect. It also requires timing the measure of drug effect so that the peak effect is not missed. For pharmacokinetic and pharmacodynamic studies using continuous infusions rather than boluses, we encourage investigators to calculate and report the time of peak effect predicted in their patients. We also encourage clinical investigators and those simulating drug concentrations to use the t^{peak}method when combining pharmacokinetic and pharmacodynamic parameters from separate studies.

The argument that C^{e}cannot have two or more local maxima was invented by D. Russell Wada, Ph.D. (Pharsight Corporation, Mountain View, California).

##### Appendix

In this Appendix, we discuss the details of the theory and implementation required to compute t^{peak}from k^{e0}and *vice versa* . As before, let the concentration of drug in the plasma following a bolus dose with no drug already in the system have the polyexponential form given by equation 1. Time t is always taken to be nonnegative, with the bolus occurring at time 0. The only assumption we make about the parameters in equation 1is that A^{1}, …, A^{n}and λ^{1}, …, λ^{n}are all positive. We can assume without loss of generality that λ^{1}, …, λ^{n}are distinct; if some of the λ^{i}are equal, we can collect terms in equation 1to obtain distinct λ^{i}(with a corresponding reduction in n). Since the λ^{i}are distinct, we can label A^{1}, …, A^{n}and λ^{1}, …, λ^{n}so λ^{1}> … > λ^{n}> 0.

Also as before, we let the concentration of drug in the effect site be related to the concentration in the plasma by a first-order process with rate constant k^{e0}, which we assume to be positive. Then the concentration in the effect site resulting from the bolus is the convolution of C^{p}(t) and k^{e0}exp(−k^{e0}t), giving equation 2. We must assume k^{e0}≠λ^{i}for all i.

Let t^{peak}be the time at which C^{e}reaches its maximum value. That is, t^{peak}is defined by C^{e}(t^{peak}) > C^{e}(t) for all t ≠ t^{peak}

We begin by proving an important prerequisite to further discussion.

###### Theorem 1

The time t^{peak}that maximizes C^{e}(t) exists and is unique.

###### Proof.

The effect site and plasma concentrations are related *via* the differential equation where the dot denotes differentiation with respect to t. Taking the derivative with respect to t on both sides gives Referring to equation 1, C^{p}is a sum of decreasing functions and is therefore decreasing, which means that Ċ^{p}(t) is negative for all t. Let t* be any critical point of C^{e}, *i.e.* , a point at which the derivative of C^{e}equals 0 (assuming, temporarily, that C^{e}has at least one critical point). Then so that t* is a local maximum of C^{e}. Since all local minima, local maxima, and saddle points of C^{e}must occur at critical points, we can conclude that C^{e}has no local minima or saddle points, only local maxima. If C^{e}has two local maxima, there must be a local minimum between them, which is a contradiction because there are no local minima. Therefore, C^{e}has at most one local maximum and no local minima or saddle points.

According to the foregoing argument, C^{e}could still have no local maxima (*i.e.* , C^{e}could have no critical points), but C^{e}is nonnegative (since it is the convolution of two nonnegative functions) and equals zero at t = 0 and t =+∞. Thus, the only way for C^{e}to have no local maxima would be for C^{e}to be zero for all t, which is impossible because of our assumptions that k^{e0}≤ 0 and A^{i}≥ 0 for all i.

To summarize, we have shown that C^{e}increases from its initial value of zero at t = 0 to its maximum value at t = t^{peak}, then decreases back to zero as t →±∞.

In the use of t^{peak}as discussed in the main body of this article, we are presented with two problems: finding t^{peak}given k^{e0}and finding k^{e0}given t^{peak}. (In both problems we take A^{1}, …, A^{n}, λ, …, λ^{n}to be given.) We consider the first problem first.

From theorem 1, we know that the function C^{e}has a unique maximum, which lies somewhere between t = 0 and t =+∞. Thus, t^{peak}can be found using any one of several standard algorithms for finding a simple extremum of a function of a single variable 17(e.g., Ch. 10 of reference 17). All such algorithms require either a starting estimate of t^{peak}or a pair of values between which t^{peak}is known to lie. In the following theorem, we obtain values that are guaranteed to bracket t^{peak}. If only a single starting estimate is needed, any number between the bracketing values can be used.

###### Theorem 2

The time t^{peak}satisfies the inequalities (log k^{e0}− log λ^{1})/(k^{e0}−λ^{1}) ≤ t^{peak}≤ (log k^{e0}− log λ^{n})/(k^{e0}−λ^{n}).

###### Proof.

From equation 2we can write where

i = 1, …, n. The function C^{e,i}is what the effect site concentration would be if C^{p}(t) were given by just the monoexponential function A^{i}e^{λit}which is a special case of equation 1. Therefore, theorem 1 applies to C^{e,i}, so we know that C^{e,i}increases from 0 at t = 0 to its unique maximum at t = t^{peak,i}and then decreases back to 0 at t =+∞.

The time t^{peak,i}at which C^{e,i}reaches its maximum is the time at which the derivative of C^{e,i}goes to 0. The derivative of C^{e,i}(t) with respect to t is MATHwhich equals 0 at t^{peak,i}= (log k^{e0}− log λ^{i})/(k^{e0}−λ^{i}).

Consider the function f(λ) = (log k^{e0}− log λ)/(k^{e0}−λ). The derivative of f is

Since for any x it holds that log x ≤ x − 1 (with equality only at x = 1), we have that f′(λ) ≤ 0 [with equality only at λ= k^{e0}, which we must exclude because f(k^{e0}) is undefined], whence f is a decreasing function of λ. Recalling that λ^{1}> … > λ^{n}, it follows that t^{peak,1}= f(λ^{1}) < … < t^{peak,n}= f(λ^{n}).

Next, we consider what happens at t^{peak,1}and t^{peak,n}. At t = t^{peak,1}, the function C^{e,1}is at its maximum, whereas C^{e,2}, …, C^{e,n}are all still increasing (because t^{peak,1}is the smallest of the t^{peak,i}). That is, Ċ^{e,1}(t^{peak,1}) = 0, C^{e,1}(t^{peak,1}) = 0, whereas Ċ^{e,1}(t^{peak,1}) < 0, for i = 2, …, n. Then by equation 4,

whence C^{e}is increasing at t = t^{peak,1}.

By a similar line of reasoning, at t = t^{peak,n}, the function C^{e,n}is at its maximum, whereas C^{e,1}, …, C^{e,n}^{−}^{1}have all begun decreasing (because t^{peak,n}is the largest of the t^{peak,i}). That is, Ċ^{e,n}(t^{peak,n}) = 0, whereas Ċ^{e,i}(t^{peak,i}) < 0 for i = 1, …, n − 1. Therefore, Ċ^{e}(t^{peak,n}) < 0, whence C^{e}is decreasing at t = t^{peak,n}.

We have shown that C^{e}is increasing at t = t^{peak,1}and decreasing at t = t^{peak,n}, which means t^{peak,1}< t^{peak}< t^{peak,n}(except that t^{peak}= t^{peak,1}if n = 1), completing the proof.

Theorem 2 concludes our discussion of finding t^{peak}given k^{e0}, which brings us to the inverse problem of finding k^{e0}given t^{peak}. From theorem 1, t^{peak}is defined by Ċ^{e}(t^{peak}) = 0. That is, differentiating equation 2, t^{peak}is defined by equation 3. Factoring out k^{e0}, the solutions to equation 3are k^{e0}= 0 (which we do not allow) and the solutions (if there are any) to The following theorem shows that equation 5has a unique solution for k^{e0}that yields a maximum effect site concentration at time equal to t^{peak}.

###### Theorem 3

Given t^{peak}(and given A^{1}, …, A^{n}, λ^{1}, …, λ^{n}), there exists a unique value of k^{e0}such that the maximum of C^{e}(t) occurs at t = t^{peak}.

###### Proof.

By the above discussion, it suffices to show that equation 5has a unique solution. Denote the left-hand side of equation 5by f(k^{e0}). Then

which is positive. In addition, k^{e0}e^{−keotpeak}tends to zero as k^{e0}→+∞, so for all sufficiently large k^{e0},

Therefore, f(k^{e0}) < 0 for all sufficiently large k^{e0}. Thus, there is at least one solution to the equation f(k^{e0}) = 0.

Whereas theorem 1 states that given k^{e0}there exists a unique t^{peak}, theorem 3 states the converse. Together, theorems 1 and 3 prove that there is a one-to-one relationship between t^{peak}and k^{e0}. Moreover, by the proof of theorem 3, t^{peak}is a decreasing function of k^{e0}, so that the inverse function, *i.e.* , k^{e0}as a function of t^{peak}, is also decreasing. These properties agree with our intuition that as the rate constant k^{e0}increases, the time of maximum effect site concentration decreases.

Solving for k^{e0}given t^{peak}means solving equation 5for k^{e0}. Thus, k^{e0}can be found using any one of several standard algorithms for finding a zero of a function of a single variable [*e.g.* , Ch. 9 of reference 17]. All such algorithms require either a starting estimate of k^{e0}or a pair of values between which k^{e0}is known to lie. We have assumed k^{e0}to be positive, so it remains to find an upper bound for k^{e0}. Such a bound is given in the following theorem. If only a single starting estimate is needed, any number between zero and the bound given by the theorem can be used.

The proof of the theorem makes use of the following inequality for the function log x.

###### Theorem 4

If 1/2 ≤α≤ 1, log x ≤ (x − 1)α for all x ≥ 1.

###### Proof.

It is convenient to treat the case α= 1 separately. If α= 1, the inequality claimed in the theorem reduces to log x ≤ x − 1, a standard result.

To begin the proof of the theorem for the remaining cases, let α satisfy 1/2 ≤α≤ 1 and define f(x) = log x − (x − 1)α. We need to show f(x) ≤ 0 for all x ≥ 1. We henceforth take x ≥ 1. From elementary calculus,

and

Subtracting these integrals,

writing the integrand over a common denominator. To show f(x) ≤ 0, it suffices to show that the integrand in equation 6 is nonpositive for all t ≥ 1. Let g(t) = 1 −αt(t − 1)α^{−1}, the numerator of the integrand. Because the denominator of the integrand is positive for t ≥ 1, it suffices to show that the numerator g(t) is nonpositive for all t ≥ 1.

To show g(t) ≤ 0 for all t ≥ 1, we show that the maximum value of g(t) for t ≥ 1 is nonpositive. The derivative of g is

which goes to 0 only at t = 1/α. The second derivative of g is

and g″(1/a) =−(1 − a)^{α−2}/a^{a−4}≤ 0 (recalling 1/2 ≤α < 1). Therefore, the point t = 1/α is a local maximum for g. The value of g at t = 1/α is g(1/α) = 1 − (1/α− 1)^{α−1}. Since lim^{t}^{→}^{1}g(t) =−∞ and lim^{t→+∞}g(t) =−∞, the point t = 1/α is the global maximum for g for t ≥ 1.

We have shown that g(1/a) = 1 − (1/a − 1)^{α−1}is the maximum value of g(t) for t ≥ 1. The inequality 1 − (1/α− 1)^{α−1}≤ 0 holds if and only if 1/α− 1 ≤ 1, which is equivalent to α≥ 1/2. Thus, for 1/2 ≤α < 1, the function g(t) is nonpositive for t ≥ 1, completing the proof.

Using the theorem, we can establish an upper bound for k^{e0}.

###### Theorem 5

The effect site rate constant k^{e0}satisfies k^{e0}≤ max(1/t^{peak},λ^{n}+ 1/(λ^{n}t^{peak}^{2})).

###### Proof.

From theorem 2, t^{peak}≤ (log k^{e0}− log λ^{n})/(k^{e0}−λ^{n}). We consider two cases, k^{e0}≤λ^{n}and k^{e0}> λ^{n}.

First suppose k^{e0}≤λ^{n}. We can write

where the inequality follows from the fact that log x ≤ x − 1 for any x. Therefore, t^{peak}≤ 1/k^{e0}, and reciprocating both sides, k^{e0}≤ 1/t^{peak}.

If instead, k^{e0}> λ^{n}, then we can write

applying the theorem with α= 1/2. Thus, we have

and solving this inequality for k^{e0}yields k^{e0}< λ^{n}+ 1/(λ^{n}t^{peak}^{2}).

We have shown that either k^{e0}≤ 1/t^{peak}or k^{e0}≤λ^{n}+ 1/(λ^{n}t^{peak}^{2}), so k^{e0}cannot exceed the greater of these two upper bounds.