Mixture Model for Individual and Combined Data Using Estimating Equations

When performing analysis of individual data on the application of a particu lar drug, it is useful to study the within variability. But when two drugs are used in combination, it is of more interest to study any combination effects on the subjects. In this paper we consider a new analytical framework that is a combination of the individual and combined data analyses, based on an estimating equation approach. The proposed analyses utilize a stochastic model for a two-drug combination and derive the mean and the variance terms based on Ito’s calculus. The proposed estimation methods are used to estimate model parameters from both individual and combined data, and they provide the basis for model free synergy tests. The strength of the fit of the model to the data is examined by statistical measures and the graphical method. Simulation studies were performed to show the strengths of the proposed approach in estimat ing the model parameters. A synergy test of the model fitted by the individual subjects confirmed that the combination of the isomers under study is synergistic in nature.


Introduction
Co mbin ing drugs in biological systems is a common practice, and the effect can be one of three types: additive, synergistic, or antagonistic. Mathematical modeling is often used to describe the relationship of two drugs. One such model, referred to as the most suitable, is the Loewe additivity model under the assumption of no interactions in the model as in [1]. Th is model refers to a co mbination of t wo concentrations of drug A and drug B being additive and equivalent to an interaction effect. Reference [2] has established an extension of a model free test for synergy in mu ltip le drug comb inations. Reference [3] used additive models generated fro m stochastic differential equations for combinations of two aesthetic agents. In this paper, we use the later approach for making the model for data involving two isomers applied in co mbination. We propose the necessary and sufficient conditions to validate the reference [2] test.
It is well-known that the flow of a chemical in the hu man body and its concentration at various times can be described by d ifferent ial equat ions . A fter an in it ial dose o f t he chemical is in jected int o th e s ystem, so me o f it will escapefrom the system and some will remain to decay over time. Observat ion o f th is pheno menon suggests that the concentration of the chemical in the system at time t Here, α is a positive constant that denotes the relative rate of elimination of the chemical fro m the body and is a generally decreasing function. The model is based on one-compart ment linear model with infusion. Kinetic models of this type have been extensively studied (see e.g. the book of [4], papers of [5], [6] and many other works).
As pointed out in [7], the function may be subject to random fluctuations fro m a variety of physical and physiological sources. This has led to the introduction of stochastic versions of the model (1), where the function is assumed to contain a white noise co mponent. In this case, equation (1) is more properly interpreted as a stochastic differential equation x t g dt hdw α = − + + (2) where w is a Wiener process and g and h are deterministic functions. Models of this type have been studied in [7], [8], [9], [10], [11] and others. An important feature of these works is the calculation of the internal variability of the system as defined by the variance of the process .
The objective of this paper is to develop a stochastic model to study the nature of the interaction of two isomers, S and R, acting in co mbination on a single individual. We assume that the combined concentration is an a priori unknown linear co mbination of the two isomers. We determine this linear co mbination using a synergy test.
The underlying methodology is as follows. The amounts of isome r S and iso mer R in the system at any given time are assumed to follow diffusion processes and of the form (2), driven by the same Wiener process w . We further assume that the combined concentration of the chemicals is given by the linear co mbination of the two processes, for some values of 1 c and 2 c to be determined nu merically. In Section 2, we give formu lae for the mean and variance of the process and identify the distribution of as Gaussian. The dataset is described in section 3, fo llowed by the description of the synergy test and the statistical methodologies in the subsequent sections. We present results in section 6 and a simu lation study comparison in section 7. In Section 8, we provide the discussion of the findings.

Distribution Theorems
We specialize the model (2). Assume that the level of concentration of the isomer in the system follows a linear stochastic differential equation of the form Here v is the initial amount of isomer, β is its rate of absorption, κ is a diffusion constant, and w is a standard Wiener process. The rates α , β , and the diffusion coefficient κ are considered to be constants. We assume that the initial concentration ) 0 ( x is zero. The following results will be needed in Sect ion 6. We now describe the distribution of a linear comb ination of two processes satisfying equations of the form (3). To this end, consider the linear stochastic differential equations In the sequel, let denote the process The process is Gaussian with a mean are two constants as used in equation (9). We refer the α 's and the β 's as the main parameters for the model, κ 's as the variance parameters, and c 's as the synergy parameters.

Description of Data
This study focuses on the concentrations of two stereoisomers, referred to as S and R, in the human bloodstream. Isomers are essentially identical chemical substances that differ only in their stereochemistry. That is, the elemental ma ke-up of the substances is exactly the same, but the 3-dimensional positioning is different. Specifically for S and R iso merism, certain co mponents of the compound are arranged in a different circu lar pattern about the center. S isomers indicate "left-hand" arrangements while R isomers indicate "right-hand" arrangements. These differences are often compared to the direction of spokes on a wheel as in [12].
When the experiment begins, a b lood sample is taken fro m each of seven volunteers (referred to as subjects). It is determined that the isomers are not present in any of the initial blood samp les. That is, the concentration at time ze ro is 0 ng/mL as in [13]. A mixture containing equal parts of both isomers is applied via eye drops to the subjects at a rate of one drop per minute for three minutes. At the five minute time ma rk and eleven intervals thereafter, blood samp les are taken fro m each subject and the concentration (in ng/mL) is recorded. We define the concentration of the single dose (either R or S) by the concentration amount of the isomer over eleven time points. We also define the concentration of the combined dose by using the sum o f the S and R concentrations over the corresponding eleven time points.

A Study for Synergy
Reference [14] proposed a synergy test of two-drug combinations that does not require modeling a response surface. They have used a Model Free Test (MFT) to establish a sufficient condition for synergy at a co mbination dose. The paper discusses the study of two or more agents in ( combination of dose ratios and describes synergy as it applies to dose equivalents of the mathemat ical equation of the form below. Here, x i is the dose (or concentration) of the ith drug given in a co mbinat ion of n drugs and x ie is the dose (or concentration) of the ith drug given individually that would produce the magnitude of the co mb ination concentration. If the value of the left-hand side of the equation (12) is less than 1, the effect is synergistic; if the effect is greater than 1, the effect is antagonistic. The equation (12) In the fixed rat io design, as in [15] of two drugs, individual compounds are combined together in amounts such that the proportion between them is constant. In other words for t wo levels of drugs we look at linear co mb inations of the form or which must either be greater than, less than, or equal to 1.
In the exp ression, are the functions of one drug in the absence of the other and is the correlat ion coefficient that depends on the ratio of the co mb ination of the two. For the synergy, we use the following tests, for isomer R, and H o : for isomer S, for a specified value o f r. A necessary condition for this hypothesis test is that the power for the specified value of r is at least 0.5. Also, since correlation is time dependent and expected to change as concentration combination changes, we establish a bound on r for significant results of the synergy tests for all subjects.

Statistical Methodology
The NLIN procedure in SAS uses the mean and variance calculated fro m the stochastic model without any assumptions about the parametric form of the distribution. The procedure uses the least squares method to fit the curve to the observations and estimate the model parameters. The process requires that the first derivatives of the equation with respect to each model parameter to be estimated. Mean Squared Error (MSE) is calculated as a result of convergence of the NLIN procedure.
The NLM IXED procedure in SAS requires that the mean and the variance expressions for the equations be supplied as the initial input with the assumption that model has a normal distribution with the stated parameters. The procedure uses ma ximu m likelihood estimates of the parameters wh ile fitting the model to the observations. It includes information about the mean and variance for the combined data, and hence we estimate the co mbination coefficients 1 c and 2 c of the drug efficacies. We use individual data for later analyses where 1 c and 2 c are considered to be known constants.
Akaike's information criterion (AIC) is a useful statistic for statistical model evaluation and has been widely accepted in some areas of statistics, eg. See [16]. It is calculated for each selected model as AIC = (n)ln (SSE/n) + 2k , where k is the number of parameters to be estimated and SSE stands for sum of squared errors. A low value for AIC indicates a better fit as described in [17]. The value of AIC is computed after the convergence of the NLMIXED procedure. The value of AIC is calculated for NLIN procedure fro m the respective MSE values.
As suggested in [18] we use the Wilco xon Ran k Test for the linear co mbination of hypotheses, that were described in the synergy test establishes a sufficient condition for rejecting the hypotheses at a .05 significance level. Since the distribution of the time data for each subject is unknown, the nonparametric min test is appropriate and widely used. We perform the power study using the WILCOX.TEST procedure in R for a simu lation size o f 5000 data sets. For the synergy hypotheses, we test seven pairs of hypotheses against one-sided alternatives, one pair per subject on each isomer type for co mbined and individual data.
For model (9), using R software we generate 20,000 data sets for each patient using the estimated parameters fro m both the NLIN and NLMIXED procedures. We supply the initial estimates of the parameters, and use the NLS function to check the convergence of the model parameters to their initial values.

Results
The initial dose amounts are considered known and used as constants throughout the analyses. We use the combined data for seven patients to estimate the parameters in equation (9). We use the equations (10) and (11) to input the mean and variance required in the procedure. NLMIXED was necessary to incorporate the subjects as blocking variables. This allo wed us to analyze the data in its entirety and estimate the synergy parameters 1 c and 2 c as given in Table 1. Th is capability is not available fro m the NLIN procedure. However, we use the individual data to estimate mean and variance parameters fo r seven subjects using both the NLIN and NLMIXED procedures. The results are presented in Tables 2 -3, using the given coefficients of 1 c and 2 c fro m Tab le 1.  The above two procedures generated parameter estimates with reasonably low standard errors of estimat ion. An AIC comparison further shed insights on the nature of the analysis done. The AIC values in Table 4 show lower nu mbers for the NLIN procedure in every subject than the NLMIXED and points towards a better fit of the data by the former procedure.
Using the parameter estimates obtained fro m the data in previous tables, we fit the estimated equations to the individual data as shown in Figure 1. The fitted curves below indicate an extremely good fit of the model to the data by both the procedures. A careful look at the fitted curves by the NLIN procedure confirms that it fo llo ws the data for the individual subjects slightly better than the curves fitted by the NLM IXED p rocedure.
The results for power study of each are given in Table 5. The lower bound for r is the ma ximu m (i.e. it is equal to 1.00). The upper bound for r is the estimate that generated a power, wh ich we accept (i.e. it is above 0.5). Using these ranges, we perform the synergy tests on the observed data.
In table 5 the sixth subject shows the maximu m strength of maintaining the power. The synergy tests for individuals based on the choice of correlation coefficient between .08 and .23 (over all data) gives p-values less than .05 for each subject when performed fo r both the isomers. The W-score for the W ilco xon Rank-Su m statistics and the p-values are shown in Tables 6 and 7. A description of the estimated values fo r the correlat ion coefficient is shown below in the isobologram in Figure 2. The theoretical isobole described as in [14] illustrates two hypothetical doses with a generic placement of the synergistic locations. The straight line represents the e-theoretical line o f additiv ity connecting and . In Figure 2, we use the same line setup and place the correlation bounds to show the synergistic power and possible location for the dose combination of the data.

Simulations and Efficiency Comparison for Parameter Estimations
A Bootstrapping Monte Carlo simulation study was conducted for individual subjects to compare the inferential performance of the NLIN and NLMIXED procedures described in section 5. We considered the parameter estimates fro m the data as the initial starting points. After 20,000 iterations, the estimates converged to assigned criteria and produced the average of each parameter value with the corresponding standard deviation. We used the known estimates of 1 c and 2 c in the simu lation to be consistent with the synergy test. In most cases in Tables 8 -9 the estimates are very close to the true parameters in Tables 2 -3. Differences shown in  indicate that some o f the variance ( ) parameters were incorrectly estimated by the data, though the main parameter estimates were quite close. There was ) , 0 ( ) , ( * : ,κ κ only one indication of a large difference that was detected by the simulation. Table 11, the NLMIXED difference table shows a large difference for two main parameter estimates for Subject 1. The NLIN procedure seems to converge faster than the NLMIXED procedure and the parameter estimates are also closer to the actual values obtained from the experimental data. We encountered problems of occurring singular values using the NLMIXED method. The variance estimates seemed not to converge well to the true values for the either procedure. However the bulk of the convergence was smooth, giving precise estimates with very small standard errors.

Conclusions
In this paper we have addressed one very important question of drug combination. Pharmacology studies deal with co mbined data that are interactive in nature. In mathematical terms an interaction is usually represented by mu ltip licat ion of terms. It is not quite clear for drug interactions whether the combined concentrations are additive, mult iplicative, or inh ibitive in nature. Co mbination of isomers in eye doses gives a reason to test for synergy to assert the researchers' effort to show an additive effect.
The proposed model uses very simp le stochastic differential equation techniques to solve for the model mean and variance in closed form. The two different statistical procedures are used here to confirm that our model has the flexib ility for use by practitioners.
The proposed model is simple in nature and uses a limited number of dose combinations. We were faced with the challenge of using the existing data. But the data collection can be done based on the recommendations in the paper [14].
While many other methods of drug assessment and model building are available, our model tests for a mathematical synergy instead of a biological synergy. A biological synergy may or may not happen in the experiment; but, if data are collected based on the combined dose levels as suggested by other researchers, a mathematical synergy can be tested as shown in this paper.