Emergence of power laws in the pharmacokinetics of paclitaxel due to competing saturable processes.

PURPOSE
This study presents the results of power law analysis applied to the pharmacokinetics of paclitaxel. Emphasis is placed on the role that the power exponent can play in the investigation and quantification of nonlinear pharmacokinetics and the elucidation of the underlying physiological processes.


METHODS
Forty-one sets of concentration-time data were inferred from 20 published clinical trial studies, and 8 sets of area-under-the-curve (AUC) and maximum concentration (Cmax) values as a function of dose were collected. Both types of data were tested for a power law relationship using least squares regression analysis.


RESULTS
Thirty-nine of the concentration-time curves were found to exhibit power law tails, and two dominant fractional exponents emerged. Short infusion times led to asymptotic tails with a single power exponent of - 1.57 +/- 0.14, while long infusion times resulted in steeper tails characterized by roughly twice the exponent. The curves following intermediate infusion times were characterized by two consecutive power laws; an initial short slope with the larger alpha value was followed by a crossover to a long-time tail characterized by the smaller Beta exponent. The AUC and Cmax parameters exhibited a power law dependence on the dose, with fractional power exponents that agreed with each other and with the exponent characterizing the shallow decline. Computer simulations revealed that a two- or three-compartment model with both saturable distribution and saturable elimination can produce the observed behaviour. Analogous linear models did not provide good fits over the range of values collected empirically. Furthermore, there is preliminary evidence that the nonlinear dose-dependence is correlated with the power law tails.


CONCLUSIONS
Assessment of data from published clinical trials suggests that power laws accurately describe the concentration-time curves and nonlinear dose-dependence of paclitaxel, and the power exponents provide new\ insights into the underlying drug mechanisms. The interplay between two saturable processes can produce a wide range of behaviour, including concentration-time curves with exponential, power law, and dual power law tails.


INTRODUCTION
One of the major challenges in clinical oncology is estimating the optimum dose and dosing time of an anticancer drug for a given patient. The concept of linearity in the body's handling of a drug is important, since the concentration of the drug as well as derived parameters scale simply with both dose and time. Nonlinearity, however, implies that the relationships are less straight-forward. In this study, we investigate new ways to assess and quantify nonlinear pharmacokinetic behaviour using power laws, with emphasis on their origins and applications to the drug paclitaxel.

The pharmacokinetics of paclitaxel
Paclitaxel (1) is derived from the bark of the European Yew tree, and it is active against many cancers, including solid tumors related to ovarian, breast, non-small cell lung, and head and neck cancers (2). Because it is poorly water-soluble, the current formulation is in a solvent of Cremophor EL (CrEL) and dehydrated alcohol. Paclitaxel is typically administered by intravenous infusion over 1, 3, 6, or 24 hours. _______________________________________ Because of patient reactions to CrEL, alternative formulations of the drug have been introduced, including Genexol-PM (3) and ABI-007 (4).
Paclitaxel is eliminated predominantly through metabolism in the liver by cytochrome P-450 enzymes (5). It works by binding to and stabilizing the microtubules within cells, especially dividing cells such as tumor cells, eventually leading to the inhibition of cell replication resulting in cell death due to apoptosis (6). Paclitaxel has a long residence time within the body and can stay trapped in cancer cells for over a week (7). Paclitaxel is also highly bound to CrEL micelles, plasma proteins, platelets, and red blood cells (8).
Clinical trials have revealed that the area under the plasma curve (AUC) and the maximum plasma concentration (C max ) for paclitaxel increase disproportionately with an increase in dose. Not surprisingly, compartmental models with linear rate constants have provided less than adequate fits to paclitaxel concentration-time curves. As a result, two-and three-compartment models with both saturable distribution and saturable elimination have been used to model the clinical data (9,10,11). Figure 1 shows the central plasma compartment (#1), a saturable binding compartment (#2), and an optional linear binding compartment (#3). The saturable distribution has been attributed to either transport (11) or binding (12) processes. The pharmacodynamic effects of paclitaxel correlate best with the duration of time that the plasma concentration remains above a critical value, estimated as 0.05 μM (10).

Power Laws in Pharmacokinetics
The variable y follows a power law function of x if where the exponent b determines the shape of the relationship between y and x, and a controls the magnitude of the gain in y. If |b|< 1, y increases or decreases more slowly than x, and if |b|> 1, y changes more rapidly than x. Equation (1) possesses the property of scaling whereby if x is multiplied by a factor L, the constant a changes but y remains proportional to x b : The behaviour of y is said to be self-similar over the range of x for which this relationship holds. Figure 1. The three-compartment model with both saturable distribution and saturable elimination from the central compartment, as well as an optional linear binding compartment.
The solution to a compartmental model with constant coefficients takes the form of a linear superposition of exponential terms, and the resulting concentration-time curve exhibits an exponentially-decaying tail. However, there is evidence that the concentration-time curves for many drugs exhibit long-time power law tails of the form where T marks the onset time and α is the power exponent. Negative power laws were first applied, empirically, to describe the washout of boneseeking radioisotopes (13,14). Subsequently, other types of clearance curves have been fit by a single power law, two sequential power laws, or the gamma function, y(t)=at -α e -βt (15,16,13,17). Although most of these applications of the power law have been strictly empirical, several models have been proposed to explain the shape of the curve, including a stochastic random walk model based on the cycling of molecules in and out of the plasma (18), a set of convection-diffusion equations for transit in the liver (13), gamma-distributed drug residence times (19), and fractal elimination kinetics (20). Certain types of kinetics cannot, by themselves, generate power law tails. Consider a one-compartment model post-infusion, where dC/dt = -kt. In the case of classical kinetics, k=k 0 , and the resulting concentration-time curve has an exponential tail. In the case of fractal kinetics, k=k 0 t -h , and the curve has a stretched exponential tail. In the case of Michaelis-Menten kinetics, dC/dt=-v max C/(K M +C), and the curve exhibits an initial linear segment followed by an exponential tail at low concentrations. In this paper, we demonstrate how a multi-compartment model with two Michaelis-Menten processes can produce both single and dual power law tails.
Power laws have also been used to model the dose-proportionality of pharmacokinetic parameters, such as the area under the curve, AUC, and the maximum concentration, C max . If a doubling of the dose of a drug produces a doubling in one or more pharmacokinetic parameters, the system is considered to be dose-proportional and linear. If the parameter decreases or increases by a factor other than two, the system is considered to be not dose-proportional and nonlinear. To quantify the extent of nonlinear dose-dependence, Gough et al. (21) proposed the "power model": where P represents the pharmacokinetic parameter and D is the dose. However, the authors only considered the cases corresponding to β=0 and β=1, and their only suggestion for fitting a nonlinear relationship was the addition of higherorder polynomial terms to Eq. (4). In addition, the authors stress that their model is empirical, and their application of the power law is not mechanistic. To the best of our knowledge, a link has yet to be made between power law concentration-time curves and power law dosedependence. Here, we report non-integer values of β and attempt to relate them to characteristics of the underlying concentration-time curve.

METHODS
A power law is best identified through a log-log plot, since taking the logarithm of both sides of Eq.
(1) yields the linear relationship log(y)=blog(x)+a, where the slope of the line is equal to the power exponent b. Thus, the existence of a power law can be tested for by performing regression analysis on log-transformed data, with the goodness of fit being evaluated using the R 2 metric (a value of 1 corresponding to a perfect fit). Forty-one sets of concentration-time data from 20 published clinical trial studies were digitized and inferred using Macromedia Fireworks Version 4. The data were tested for power law tails of the form expressed by Eq. (3), where T coincides with the end of the infusion. The mean number of data points in the tail was 9±2 (minimum of 6 and maximum of 13). Eight sets of AUC and C max data were taken directly from 10 published studies and were fit to Eq. (4). The mean number of dose levels per study was 4±1 (minimum of 3 and maximum of 6). Appendix A is provided for the ease of reference to the reader and contains a collection of data tables summarizing the key pharmacokinetic data in the publications analyzed by the present authors. Figure 1 shows the central plasma compartment, transfer to a second compartment via Michaelis-Menten kinetics, and optional transfer to a third compartment via linear kinetics. Saturable distribution has been attributed to either transport (11) or binding (12) processes. The diffusion equations describing the model are given by: and . 3 31 The superscript "d" indicates parameters that describe the distribution process, the "e" superscript indicates the elimination process, Simulations of the model illustrated in Fig.  1 and described by equations (6) to (8) were performed using code written in C++ to numerically solve the set of ordinary differential equations using a fourth-order Runge-Kutta algorithm (22). Data points corresponding to a specific set of parameters were then generated, and the AUC was calculated using Simpson's method.

Power law tails
Thirty-nine of the 41 concentration-time curves exhibited power law tails. The values calculated for the power exponent α are listed in Table 1.
Although the exponent was relatively independent of the patient characteristics (such as weight, age, sex, and the type and stage of cancer) and the dose level, it varied markedly with the length of the infusion. For short infusions (1-h duration), a single long-time tail was observed with a power exponent of α =1.57±0.14. The tails persisted up to 24 h and in one case up to 36 h. For long infusions (6-h or 24-h duration), a single long-time tail was also observed but with an exponent of over 3. The tails extended up to 24 h for the 6-h infusions and up to 48 h for the 24-h infusions.
After confirming the existence of power law tails in paclitaxel concentration-time curves, we wanted to test whether they could be generated by the competition between two saturable processes. First, we confirmed that the data in table 1 was best fitted by the model illustrated in figure 1. We indeed found that the Akaiki Information Criterion (AIC) (23) was lower for both the two-compartment version (AIC = -16.3 ± 3.6) and the threecompartment version (AIC = -15.4 ± 2.9) than for either a linear two-compartment model (AIC = -5.8 ± 1.8) or a two-compartment model with fractal Michaelis-Menten elimination kinetics (24) (AIC = -6.5 ± 2.1).  (31,11,37).
In contrast, for the intermediate infusion time of 3 h, the concentration-time curves exhibited a crossover between two power law regimes. At the end of the infusion, there was an initial 1 -2 h slope whose power exponent was equal within error to that exhibited by the 6-h and 24-h infusion curves. This initial slope was followed by a longtime tail with a power exponent equal within error to that exhibited by the 1-h infusion curves. The length of the long-time tail ranged from 9 to 68 h post-infusion. Figure 2 shows a plot of 3-h infusions of paclitaxel, in which the dual power law nature of the tail is evident and independent of dose. Figure 2(a) shows the semi-log plot of the data illustrating that it cannot be fitted to an exponential curve. Figure 2(b) shows clearly that the data is linear in the log-log plot suggesting an asymptotic power law decay. Figure 2(c) illustrates that a linear version of equations 6 to 8 does not fit the data and the nonlinear saturable compartments are required. The rise of each curve also appears to follow a power law relationship. These curves demonstrate a flip-flop situation because the initial slope is steeper than the terminal long-time tail. One would expect the opposite, since smaller plasma concentrations should be cleared more rapidly. Usually, a flip-flop is linked to orally-administered drugs whose absorption is the rate-limiting step (38). In the case of paclitaxel, however, the flip-flop is associated with the tissue distribution process. The observed power law behaviour is likely not a consequence of binding to the CrEL formulation vehicle or blood components. Analysis of two CrEL-free formulations of paclitaxel produced results consistent with those observed for regular paclitaxel. We conjecture that the steep curves correspond to the situation where the distribution process is not saturated, allowing the maximum fraction of drug to be distributed outside the plasma. The power exponent, therefore, reflects the maximum temporary and permanent transfer of drug outside the plasma. On the other hand, the shallow curves result when the distribution process is saturated, and the power exponent predominantly reflects the elimination process. This occurs when the drug is infused too quickly. The fact that the power law tails persist even at low concentrations (below the K M values) provides additional information about the system. This continued adherence to a power law indicates a failure of the drug to attain a steady state in the peripheral compartment(s). Therefore, drug distribution to and release from the tissues plays a dominant role in the pharmacokinetics of paclitaxel at all plasma concentrations.

Simulations
Since the linear model was unable to adequately approximate the experimental data (see Fig. 2c), it is clear that a non-linear model is required. After confirming the existence of power law tails in paclitaxel concentration-time curves, we intended to test whether they could be generated by the competition between two saturable processes. We performed simulations of the model shown in Fig. 1 using the parameter values reported by Kearns et al. (9) (summarized in Table 2) and found that the model can indeed reproduce the observed behavior. Figure 3 shows the shape of the concentration-time curve as a function of the infusion time. Experimentally, the dual power-law curve is seen in only the three-hour infusion, but the simulations shown in Figure 3(a) show a dual power-law decay for the 1-h, 3-h, and 6-h infusions. However, in the case of the 11-h infusion, the slopes differ very little within experimental error. In the 6-h case, the transition range is very short and could easily be missed in the experimental measurements. Figure  3(b) illustrates that the linear model cannot simultaneously represent the infusion and elimination parts of the curve properly. The reader should note the different scales on the time axes in Figures 3(a) and (b). The linear solutions approach zero more rapidly than the non-linear solution to the concentration equation.
For the 3-h infusion, a dual power law is evident. As the infusion time decreases, the initial slope becomes less steep, and a single power law may be observable after short infusion times under certain conditions. Conversely, as the infusion time increases, the initial slope becomes shorter and essentially unobservable, resulting in a single steep power law for the 24-h infusion curve.
An advantage of computer simulations is that they allow us to study the effect of different parameters on the shape of the curve.

b)
The duration of the initial slope is determined by the parameters describing the saturable distribution process, v d max and K d M , while the value of the initial exponent is determined predominantly by the parameters describing the saturable elimination process, v e max and K e M . The slope of the terminal curve is mainly determined by v d max , where an increase in its value produces an increase in α. The inclusion of the third, linear binding compartment was found to have an important effect. For low values of k 13 (weak binding), the slopes of the two segments remain unchanged, but the duration of the initial slope increases (Fig. 4). As k 13 increases, however, the long-time tail eventually becomes exponential instead of power law. Therefore, an increase in the strength of the linear binding process decreases the plasma concentration and minimizes the importance of the saturable processes.
In the case of paclitaxel, K d M << K e M , and therefore the distribution process saturates before the elimination process. But what occurs when the reverse, K d M >> K e M , is the case? An example is shown in Figure 5 for the two-compartment model, and a much wider range of behavior can result. Simply changing the value of v d max can produce a dual tail with a terminal exponential, a flip-flop situation with a dual power law tail, a single power law tail, a dual power law tail with a steep terminal portion, and a dual tail with a steep long-time exponential. Moreover, in Figure 5(b), we have illustrated the failure of the linear model to exhibit the flip-flop property. A similar transition between different regimes also occurs when the dose increases, the volume of distribution decreases, or k 21 increases. Therefore, this situation is much more sensitive to changes in the dosing regimen or patient characteristics. Although this discussion is not relevant to paclitaxel, the model employed to describe its behavior is more general and could be extended to other chemotherapy drugs or other formulations of taxol.

b) a) Power law dose-dependence
The dose-dependence of the AUC and C max were found to be non-linear, in agreement with the literature. The relationship in Eq. (4) provided a good fit to the data, and the results are listed in Table 3. The two values of β agree within error with each other and with the exponent characterizing the shallow long-time tails, although the latter might be a coincidence. Note that these results are valid over the therapeutic dose ranges considered, and caution should be used in extrapolating beyond this range.
Additional simulations were performed to investigate whether a relationship exists between β and α, such that the non-linear dose-dependence of the AUC and C max is correlated with the nonexponential shape of the long-time tail. Figure 6 shows that β increases as α decreases; therefore the dose-dependence becomes increasingly disproportionate as the long-time tail becomes less steep. In other words, the greater the contribution that the tail portion of the curve makes to the overall area under the curve, the greater the increase in the AUC with an increase in dose.  (10,9,24,25,29,37,32,35). b From (10,9,24,25,29,1,35,38).

DISCUSSION
We have shown that two competing saturable processes can generate concentration-time curves with power law tails of different forms. An advantage of this proposed mechanism is that it does not place a restriction on the value of α, in contrast to models suggested for the bone-seeking elements, where -1 < α < 0 (13), and for fractal kinetics, in which α=1-d s /2 and d s < 2 (40).
Before embarking on a full pharmacokinetic analysis, the tails of concentrationtime curves can be checked for a power law fit. The existence of a power law tail, especially a shallow one, can signal dose-dependent behaviour. In the case of paclitaxel, a steeper decline is more desirable, since it implies a greater distribution to the tissues as well as a decrease in the extent of non-linearity. This model is also consistent with the observation that the pharmacodynamic effect does not correlate with the AUC or C max , since the plasma concentration is not necessarily indicative of the amount of drug still present in the body.
Simulations can also be used to investigate the dose-dependence of AUC and C max beyond the current clinically-relevant range. Both parameters exhibit three distinct regions of dose-dependence on a log-log plot (Figure 7). In both cases, the initial and terminal regions are characterized by approximately β = 1 (Table 4), indicating linear kinetics when the system is well below or well above the saturable concentration range. The intermediate regions, however, are characterized by a transition to a nonlinear regime with β > 1. The curves in Figure 5 fall within this dose range. Because the onset dose of the intermediate slope is higher for AUC than for C max , and the slope persists over a longer range, there appears to be a lag between the occurrence of disproportionately higher maximum concentrations and an overall noticeable effect on the shape of the concentration-time curve. In summary, power law analysis is helpful in predicting non-linear kinetics, and Figure 7 emphasizes that the concept of linearity is valid only over a specified dose range.
Finally, power laws may indicate self-similarity, as discussed above. This may pertain to both temporal and spatial self-similarity. Power law relationships in concentration-time curves may imply that pharmacological processes are linked over different time and/or size scales. To date, allometric scaling has been applied to pharmacokinetics to both assess inter-species and intra-species variation in pharmacokinetic parameters. However, the scaling of pharmacological processes within an individual over different time or length scales has yet to be discussed and warrants further investigation. Of particular interest in the role of fractal architecture of organs such as the liver or kidney and the vascular networks in giving rise to power-law behavior in pharmacokinetic and metabolic processes.

CONCLUSIONS
While an examination of the exact mathematics underlying a model with two saturable processes is beyond the scope of this study, it will be performed in depth in a follow-up article. Here, we wish to bring attention to the existence of the power laws and the value of the power exponents. To the best of our knowledge, this is the first study to report power law tails in the concentration-time curves of an anticancer drug, as well as the first study to relate the existence of power law tails to saturable processes. Although a single saturable process cannot produce a power law tail, two competing saturable processes can produce a range of behaviour that includes single and dual power law tails. Furthermore, although power law tails have been reported for clearance curves, this is the first study to investigate the tails of infusion curves and their dependence on the infusion rate. The exponent of the power law tail can provide insight into the underlying drug processes. For paclitaxel, only two exponents were found to characterize the behaviour following short, intermediate, and long infusion times. Analysis of the power law tails can provide valuable information for clinical applications. For example, longer infusions leading to a steeper tail are actually best for paclitaxel, because they indicate that a greater portion of the drug is being transferred to the tissues.