#### ArticlePlus

*e.g.*, propofol, isoflurane) and an analgesic (

*e.g.*, fentanyl, nitrous oxide). Anesthetic drugs are often combined because they interact synergistically to create the anesthetized state.

^{1}The combination of two doses (d

_{1}and d

_{2}) can be represented by a point on a graph, the axes of which are the dose axes of the individual drugs (fig. 1). The isobole connects isoeffective doses of the two drugs when administered alone, D

_{1}and D

_{2}. If the isobole is straight (fig. 1A), then the relation is additive. If the isobole bows toward the origin (fig. 1B), then smaller amounts of both drugs are needed to produce the drug effect when administered together, so the relation is supraadditive or synergistic. If the isobole bows away from the origin (fig. 1C), then greater amounts of both drugs are needed to produce the drug effect when administered together, so the relation is infraadditive.

^{2}Response-surface methodology is used for two principal purposes; to provide a description of the response pattern in the region of the observations studied and to assist in finding the region in which the optimal response occurs. Our model is a straightforward extension of the sigmoidal concentration–response relation for individual drugs. We test the proposed model using data from a study of the interaction of midazolam, propofol, and alfentanil with loss of consciousness.

^{3}

#### Methods

##### Model Development and Mathematics

_{0}is the baseline effect when no drug is present, E

_{max}is the peak drug effect, C

_{50}is the concentration associated with 50% drug effect, and γ is a “sigmoidicity factor” that determines the steepness of the relation.

*e.g.*, dose, plasma concentration, or area under the curve). For models of probability, such as the probability of moving in response to surgical incision, E

_{0}is 0 and E

_{max}is the maximal probability (usually assumed to be 1).

_{50}

^{γ}, we obtain an alternate form:MATH 2

*i.e.*, B/(A + B), called θ herein) of the two drugs as behaving as a new drug. This new drug, which is actually a fixed ratio of the two drugs, has its own sigmoidal concentration–response relation, as shown in figure 3. This is the basic premise of our interaction model. The mathematics are simply an extension of the model for a single drug to a model that considers each ratio of two drugs as a drug in its own right.

Equation 3 Image Tools |
Equation 4 Image Tools |

_{50}, and express the results in units (U) of potency. MATH 3 where U

_{A}is the normalized concentration of drug A, and U

_{B}is the normalized concentration of drug B. We can define a family of “drugs,” each being a unique ratio of U

_{A}and U

_{B}. Each drug will be defined in terms of θ, where θ is defined as MATH 4

_{A}+ U

_{B}. We can extend equation 2 to describe the concentration–response relation for any ratio, θ, of the two drugs in combination:MATH 5 where θ is the ratio of the two drugs, the drug concentration is U

_{A}+ U

_{B}, γ(θ) is the steepness of the concentration–response relation at ratio θ, U

_{50}(θ) is the number of units (U) associated with 50% of maximum effect at ratio θ, and E

_{max}(θ) is the maximum possible drug effect at ratio θ. Because E

_{max}, C

_{50}, and γ in equation 2 have been replaced by functions of θ, each ratio has the potential to have its own E

_{max}, C

_{50}, and γ. This allows each ratio of drug A and drug B to behave as its own drug, with its own sigmoidal concentration–response relation, which is the basic premise of the model.

_{50}(θ)” is the potency of the drug combination at ratio θ relative to the normalized potency of each drug by itself. This requires careful explanation. Let us assume that only drug A is present, in a concentration of C

_{50,A}. In this case, the drug effect is half of the maximal effect, U

_{A}= 1, U

_{B}= 0, θ = 0, and the drug concentration is U

_{A}+ U

_{B}= 1. Because we have 50% of the maximum drug effect, and 1 unit of drug, then the number of units associated with 50% drug effect when only drug A is present, U

_{50}(0), must be 1. Similarly, let us assume that only drug B is present and the concentration of drug B is C

_{50,B}. In this case, the drug effect is half of the maximal effect, U

_{A}= 0, U

_{B}= 1, θ = 1, and the drug concentration is U

_{A}+ U

_{B}= 1. Because we have 50% of the maximum drug effect, and 1 unit of drug, then the number of units associated with 50% drug effect when only drug B is present, U

_{50}(1), must again be 1. By definition, if only drug A or drug B is present, U

_{50}(θ) = 1.

_{A}= 0.5, U

_{B}= 0.5, θ = 0.5, and the drug concentration is U

_{A}+ U

_{B}= 1. If this causes 50% of maximum effect, then the drugs are simply additive at θ = 0.5, and U

_{50}(0.5) = 1. However, if this combination produces more than a half-maximal effect, then 1 unit of this combination, at θ = 0.5, is more potent than 1 unit of either drug alone (

*i.e.*, synergistic). In this case, U

_{50}(0.5) < 1. Conversely, if this combination produces less than a half-maximal effect, then 1 unit of this combination, at θ = 0.5, is less potent than either drug alone (

*i.e.*, infraadditive). In this case, U

_{50}(0.5) > 1. Thus, U

_{50}(θ) is the potency of the combination compared with the potency of either drug alone, which is 1 by definition.

_{50}(θ) are not concentration units, but rather the number of units, at ratio θ, associated with 50% of maximal drug effect. U

_{50}(θ) is 1 for θ = 0 and θ = 1. For all values of θ between 0 and 1 (

*i.e.*, all possible ratios of the two drugs), U

_{50}(θ) assumes a value determined by the data. If this value is 1, then the interaction is additive at θ. If the value is less than 1, then the drug effect is synergistic at θ. If the value is greater than 1, then the interaction is antagonistic at θ.

_{max}(θ), U

_{50}(θ), and γ(θ). Our choice is to use functions that are capable of taking a variety of shapes, so that good approximations to the true relations can be determined empirically. To provide these flexible functions we chose fourth-order polynomials of the form MATH 6 where f(θ) is E

_{max}(θ), U

_{50}(θ), or γ(θ). The coefficients (β

_{0}, β

_{1}, β

_{2}, β

_{3}, β

_{4}) are model parameters that are either constrained by the model or estimated from the data. Fortunately, two of these terms, β

_{0}and β

_{1}, can be replaced by other terms already defined.

_{max}(θ), U

_{50}(θ), and γ(θ) when only drug A is present, E

_{max,A}, U

_{50,A}, and γ

_{A}, respectively. Note in equation 6 that when θ = 0 (only drug A is present), f(0) = β

_{0}. Therefore, when f(θ) is E

_{max}(θ), U

_{50}(θ), or γ(θ), β

_{0}must be E

_{max,A}, U

_{50,A}, and γ

_{A}, respectively.

_{max}(θ), U

_{50}(θ), and γ(θ) when only drug B is present, E

_{max,B}, U

_{50,B}, and γ

_{B}, respectively. Referring again to equation 6, when θ = 1 (only drug B is present), f(1) = β

_{0}+ β

_{1}+ β

_{2}+ β

_{3}+ β

_{4}. We can rearrange this as β

_{1}= f(1) − β

_{0}− β

_{2}− β

_{3}− β

_{4}. Thus, when f(θ) is E

_{max}(θ), U

_{50}(θ), or γ(θ), β

_{1}must be E

_{max,B}− E

_{max,A}− β

_{2,Emax}− β

_{3,Emax}− β

_{4,Emax}, U

_{50,B}− U

_{50,A}− β

_{2, U50}− β

_{3,U50}− β

_{4,U50}, or γ

_{B}− γ

_{A}− β

_{2,γ}− β

_{4,γ}, respectively.

_{max}(θ), U

_{50}(θ), and γ(θ) as functions of θ. The equation for E

_{max}(θ), using the substitutions previously mentioned for β

_{0}and β

_{1}, is MATH 7

_{50,A}and U

_{50,B}, [equivalent to U

_{50}(θ) and U

_{50}(1)], are both 1 by definition. Thus, when f(θ) = U

_{50}(θ), the values of β

_{0}and β

_{1}in equation 6 are 1 and −β

_{2}−β

_{3}−β

_{4}, respectively. Therefore, the equation for potency as a function of θ can be simplified to MATH 8

Equation 10 Image Tools |
Equation 11 Image Tools |

_{2,U50}is 0, then the value of U

_{50}(θ) will be 1 for all values of θ. This means that the interaction will be additive. If β

_{2,U50}is a positive number, then U

_{50}(θ) will be less than 1 for all values of θ between 0 and 1. The effect is to magnify the term MATH 10 in equation 5, making it appear that there is more drug present. This will produce a greater than additive effect,

*i.e.*, synergy. If β

_{2,U50}is a negative number, then U

_{50}(θ) will be greater than 1 for all values of θ between 0 and 1. This reduces the term MATH 11 in equation 5, making it appear that there is less drug present. This will produce a less than additive effect. This assumes that drugs A and B have the same maximal effect. It is possible for some approaches to synergy analysis to show apparent synergy if the maximal effects of drugs A and B are not identical, even if U

_{50}(θ) = 1 for all values of θ.

_{0}and β

_{1}. The resulting equation is MATH 12

*i.e.*, β

_{2}, β

_{3}, β

_{4}) are 0. They are the equations for parabolas if the respective β

_{2}coefficient is nonzero, and β

_{3}and β

_{4}are 0. More complex shapes are generated when β

_{3}and β

_{4}are nonzero.

_{max}, U

_{50}, and γ as functions of θ for the synergistic interaction seen in figure 4. E

_{max}and γ are constant, and thus have no interaction. U

_{50}is necessarily 1 at θ = 0 and θ = 1, but is less than one between these extremes. This increases the potency of the drugs when administered in combination, resulting in the synergy seen in figure 4.

_{A}, θ

_{B}, and θ

_{C}, where MATH 13

_{A}+θ

_{B}+θ

_{C}=1. For our purposes here, we will use θ

_{B}and θ

_{C}. We again assume that for any fixed value of θ

_{B}and θ

_{C}, there is a sigmoidal relation between concentration and response. Therefore, if the three drugs could be administered to the effect site in an exactly fixed proportion, they would show a sigmoidal total concentration–response relation, where the “concentration” was the sum of the three normalized concentrations. This is precisely the notion that underlies the two-drug model. The equation for the three-drug model follows:MATH 14

_{max}, γ, and U

_{50}, are functions of θ

_{B}and θ

_{C}. The functions E

_{max}(θ

_{B},θ

_{C}), U

_{50}(θ

_{B},θ

_{C}), and γ(θ

_{B},θ

_{C}) are described in Appendix 2 (which can be found on the Anesthesiology Web site at http://www.anesthesiology.com). The important point is that instead of describing a curve, as in equations 6–10, when three drugs are present, the parameters of the sigmoidal relation (equation 12) are themselves surfaces. Figure 6 shows representative surfaces for U

_{50}as functions of θ

_{B}and θ

_{C}.

^{3}This code is available on request to the authors. ** Graphs were prepared using Excel and Mathematica (Mathematica for Windows, Version 4.0; Wolfram Research, Champaign, IL).

##### Patients and Clinical Trial Design

*et al.*

^{4}These data are also available

*via*the Anesthesiology Web site (http://www.anesthesiology.org). Dose–response relations were established after intravenous bolus doses of midazolam, propofol, and alfentanil administered individually and in combination in 400 women patients undergoing elective gynecologic surgery. Patients were assessed for hypnosis (defined as failure to open the eyes to verbal command) 4 min after midazolam administration and 2 min after propofol or alfentanil injection (the approximate times to peak effect after an intravenous bolus). When the combination being tested included midazolam, it was administered 2 min before the other drugs. The doses of midazolam, propofol, and alfentanil used and the proportion of patients achieving hypnosis for each dose category are shown in table 2. This data set was selected because it provided three two-drug combinations that could be used to demonstrate response–surface relations. In addition, the data are notable for the relatively large number of subjects and the uniformity of the experimental conditions. Lastly, it provides the ability to include a three-drug interaction model, showing the flexibility of the proposed response-surface model.

##### Statistical Analysis

_{0}was 0 by definition. In addition, we assumed that each drug was capable of causing complete hypnosis if administered in a high enough concentration.

Equation 15 Image Tools |
Equation 16 Image Tools |

_{0}(no hypnosis) and E

_{max}(θ

_{B}, θ

_{C}) (complete hypnosis) can be constrained to 0 and 1 respectively, with no interaction across the E

_{max}(θ

_{B}, θ

_{C}) parameter surface. Thus, the probability of hypnosis for any combination is MATH 15 where U

_{A}, U

_{B}, and U

_{C}are the doses of midazolam, propofol, and alfentanil, respectively. The units are fractions of each drug’s dose necessary to cause hypnosis in 50% of the population (D

_{50}). Four separate analyses based on equation 13 were performed: (1) the data for the single and paired combination of midazolam and propofol were modeled; (2) the data for the single and paired combination of midazolam and alfentanil were modeled; (3) the data for the single and paired combination of propofol and alfentanil were modeled; and (4) the complete data set for the single, paired, and three-drug combinations (table 2) were modeled simultaneously. Model parameters and standard errors were estimated using NONMEM, by maximizing the log likelihood (LL) for all 400 observations:MATH 16 where R

_{i}is the observed response of the ith individual, either 0 (respond to voice) or 1 (not respond to voice), and

*P*is the probability of response to voice for each dose combination. Equation 14 can be expressed in words as the sum of the natural logarithm of the probabilities of “no response in the nonresponders” and “response in the responders.” The contribution of the polynomial coefficients to the model was evaluated by excluding the coefficients one at a time, by determining whether the model deteriorated significantly (likelihood ratio test), and by assessment of residuals (observed

*vs.*predicted probability of hypnosis for each dose combination). Parameter estimates and standard errors were tabulated. The response surfaces for the paired interactions and the U

_{50}(θ

_{B,}θ

_{C}) parameter surface were graphed.

##### Computer Simulations

*et al.*,

^{5}Schnider

*et al.*,

^{6,7}and Scott and Stanski

^{8}were used to calculate the C

_{50}for hypnosis for midazolam, propofol, and alfentanil, respectively. This was calculated as the predicted effect-site concentration at the time of assessment after the respective D

_{50}bolus. Equation 13 was then used to calculate the predicted time course of effect (probability of hypnosis in the population) after these maximally synergistic doses of midazolam, propofol, and alfentanil, administered alone and in combination. The pharmacokinetic parameters, the calculated C

_{50}values and the time from hypnosis in 95% of the population to “no hypnosis” in 95% of the population were tabulated. The predicted time course of effect was graphed. All simulations were performed using PKPD Tools for Excel (PKPD Tools for Excel; written by C. Minto and T. Schnider, Stanford University, Palo Alto, CA).

##### Response Surfaces

#### Results

##### Midazolam–Propofol–Alfentanil Interaction

Table 3 Image Tools |
Table 4 Image Tools |
Fig. 7 Image Tools |

*et al.*,

^{4}which necessitated the exclusion of the data for 120 patients. Parameter estimates and standard errors for the analysis of the three paired drug interactions are shown in table 3. Parameter estimates and standard errors for the simultaneous analysis of the entire data set are shown in table 4. There were no significant differences in the slope parameter of the three drugs, nor were there significant drug interactions across the γ(θ

_{B}, θ

_{C}) surface. The D

_{50}values for each drug tended to be the same, regardless of the combination being modeled (tables 3 and 4). This is an expected result. Because the data for single administration was included in the interaction models, the D

_{50}values were almost entirely determined from subjects receiving each drug alone. Figure 7 shows the response surfaces for each of the paired interactions. All paired drug interactions showed significant synergy. The triple synergy parameter in the three-drug model was not statistically significant (see Appendix 2 on the Anesthesiology Web site at http://www.anesthesiology.com). This indicates that, when all three drugs are present, there is not synergy beyond that expected from the paired interactions of all three drugs.

_{50}(θ

_{B}, θ

_{C}) surface is shown in figure 8. The maximum decrease in U

_{50}(θ

_{B}, θ

_{C}) values for the paired combinations are represented by the midpoint of the three edges of the triangular surface. Values follow (expressed as a percentage decrease): midazolam–propofol = 35%; midazolam–alfentanil = 44%; propofol–alfentanil = 16%. The maximum decrease (45%) in U

_{50}(θ

_{B}, θ

_{C}) for the three-drug combination is found at the lowest point of the U

_{50}(θ

_{B}, θ

_{C}) surface, which is located at θ

_{B}= 0.144 and θ

_{C}= 0.396. This corresponds in figure 8 to θ

_{Mid–Prop}= 0.24, θ

_{Mid–Alf}= 0.47, and θ

_{Alf–Prop}= 0.73, where the point appears as a dot on the axis plane.

##### Computer Simulations

Table 5 Image Tools |
Fig. 9 Image Tools |

##### Response Surfaces

#### Discussion

_{max}equations (equation 1). Our proposed model connects these two sigmoid curves

*via*polynomial functions of θ (equation 4). Drug interactions are then sought as coefficients of the polynomials that connect the “one-drug” sigmoids to each other. Although each value of θ can have its own E

_{max}, C

_{50}, and γ, the model assumes that the underlying sigmoidal shape is preserved for all values of θ. This concept of a fixed ratio (θ) of two drugs having their own properties is not new. Carter

*et al.*

^{9}considered each drug ratio to be essentially similar to a single agent, whose response could be described by a single two-dimensional concentration–response curve. Short

*et al.*

^{4}proposed that when using combinations of sedatives “the combination should be regarded as a new ‘drug’ with individual properties, rather than merely reflecting the known properties of the individual agents.”

_{max}models, the polynomials that underlie U

_{50}(θ), E

_{max}(θ), and γ(θ) are more complex. The use of polynomial response functions to approximate complex response surfaces is common in many experimental situations. However, polynomial models with the independent variable present at powers higher than the second power are not often used because it becomes difficult to interpret the coefficients.

^{2}We proposed a very flexible model. Potentially, poor trial design, incomplete sampling of the response surface, or careless application of the curve-fitting procedures could result in estimation of parameters that provide a “good fit,” but when evaluated more closely generate an unrealistic shape for the surface. These concerns are not unique to this model and can be alleviated by proper model-selection techniques.

*e.g.*, equations 7–10) to ensure that the pharmacodynamic parameters do not extrapolate to unreasonable values. For instance, it must be ensured that U

_{50}(θ) and γ(θ) are always positive in the range of 0 ≤ θ ≤ 1. In the case of three drugs, the parameter surface should be graphed as shown in figure 8.

^{1}developed a mathematical proof showing that zero-interaction holds if and only if the isobolograms are straight lines. For this to be true, E

_{max}must be constant with respect to θ. The case of γ(θ) is more complicated. Ifγ is constant with respect to θ, there is zero-interaction if E

_{max}(θ) and U

_{50}(θ) are also constant. But ifγ

_{A}is not equal to γ

_{B}, there is no easy way to use our interaction model to test for zero-interaction (as defined by Berenbaum

^{1}). Although we agree with Berenbaum’s

^{1}treatment of the zero-interactive surface, the description of an interaction as zero-interactive, synergistic, or antagonistic may be too simplistic. For example, a drug combination can be synergistic in certain regions and antagonistic in others.

^{10}To our thinking, the focus on whether drug interactions can be reduced to simple descriptors such as “additive,” “synergistic,” or “antagonistic” completely misses the point. The interaction has the potential to be highly dimensional and complex. Rather than worry about which descriptor applies to the relation, the goal should be to characterize the response surface. From the surface one can identify the best combination to produce the specified therapeutic effect.

^{11}although we are not aware of any reason that our approach could not be combined with the more kinetic approach of indirect interactions.

_{max}model (equation 1). For example, the model could be a linear, polynomial, or exponential response function. The model could also be a bimodal response function, as seen with some hypnotics.

^{6}In fact, the model can be any function, so long as it has parameters that permit connection among the individual single-drug models. Thus, the only constraint is that the interaction model reduces to the model for drug A when θ = 0 and to the model for drug B when θ = 1.

*et al.*

^{12}These authors discuss the use of flexible nonparametric functions (splines) that are forced to obey certain constraints (for example, a spline can be constrained to resemble an E

_{max}model). Although they illustrate their approach in the context of antagonistic interactions, it can be extended to model additive and synergistic drug interactions.

^{13}