Adaptive Bayesian c-optimal design for estimating multiple ED p s under the 4-paramter logistic model

Abstract

Under nonlinear models, optimal design truly depends on the pre-specified values of model parameter. If the nominal values of the parameter are not close to the true values, optimal designs become far from optimum. In this study, we focus on constructing an optimal design that works well for estimating multiple $\text{ED}_{p}$ s taking into account the parameter uncertainty. To address the parameter dependency, an adaptive design technique is applied by incorporating Bayesian paradigm. One challenging task for the Bayesian approach is that it requires heavy computation when search for the Bayesian optimal design. To overcome this problem, a clustering method can be employed.

We propose an adaptive Bayesian c-optimal design that works fairly well for estimating multiple values of $\text{ED}_{p}$ simultaneously, whilst accounting for the parameter uncertainty. We use the flexible 4-paramter logistic model to illustrate the methodology but our approach can be extended to other types of nonlinear models. We also examine the performance of our proposed design by comparing with other traditional designs through different scenarios of simulations given a wide range of mis-specified model parameters.

Keywords

Bayesian paradigm clustering method c-optimality dose-finding study parameter uncertainty target doses estimation

1. Introduction

Experimental designs are increasingly expensive to conduct especially in dose-finding studies. A dose-finding study plays a critical role in drug development process. It mainly determines the dose of an agent for further study in large scale clinical trials. One challenge in dose-finding study is the relatively small sample size for ethical and financial reasons. Thus, setting up appropriate designs that allow researchers to gather as much accurate information as possible at minimal cost becomes crucial. One effective way to save cost and resource is to answer multiple research questions simultaneously from one single study by applying optimal design theory. However, an optimal design that works well for one research objective could perform poorly for another in dose-finding trials. Several pioneering studies have been conducted regarding multi-objective experiments recently, examples see (Miller et al., 2007; Padmanabhan & Dragalin, 2010; Hyun & Wong, 2015); While there is little in literature on developing one experimental design to investigate multiple target doses such as $\text{ED}_{p}$ in dose-finding study. The $\text{ED}_{p}$ is the effective dose that achieves $p$ % of the maximum treatment effect within the observed dose range and $p$ is given between 0 and 100 (Ting, 2006; Bretz et al., 2010). For instance, $\text{ED}_{50}$ is the median effective dose producing 50% of the maximum treatment effect. In practice, rather than studying one specific $\text{ED}_{p}$ , researchers are often interested in dual target doses such as $\text{ED}_{10}$ and $\text{ED}_{50}$ or multiple target doses such as $\text{ED}_{40}$ , $\text{ED}_{50}$ , ${ED}_{60}$ , and $\text{ED}_{80}$ to conduct the dose-range studies (Beloeil et al., 2007; Louie et al., 2001).

In general, optimal designs use different optimality criteria for studying different objectives in dose-finding trials. One primary objective for the dose-finding trial is to study an interesting target dose $\text{ED}_{p}$ . Zhang and Hyun (2016) stated that under the 4-parameter logistic (4PL) model, the c-optimal design for estimating one specific $\text{ED}_{p}$ works poorly when estimating other $\text{ED}_{p}s$ . They proposed a robust c-optimal design that works well for estimating multiple $\text{ED}_{p}s$ under the assumption that the parameter values are known. However, in practical terms, model parameters are often unknown and optimal designs under nonlinear models are very sensitive to model parameters. Typically, the construction of the optimal design is based on the initial guess of the parameter values. However, if the initial guess is not close to the true values, the optimal design is far from optimum (Chernoff, 1953; Wang & Yang, 2014). Our aim is to construct one optimal design for estimating various values of $\text{ED}_{p}$ under mis-specification in nominal values of model parameter for a biopharmaceutical application, and study the design effectiveness and robustness.

To reduce the uncertainty of the parameter values, we study an adaptive two-stage design technique by incorporating Bayesian approach into the second stage design. The Bayesian optimal design utilizes the probability distribution of the unknown parameters instead of a single set of parameters to derive a better design (Dette, 1996; Albert, 2009). A two-stage procedure offers an efficient way to learn about the parameters from the first stage and apply the accrued data to determine the augmented design at the second stage. At the first stage, a small proportion of the subjects are assigned according to a fixed design. At the second stage, the posterior distribution of the parameters is generated based on the information learned from the first stage.

One challenge for seeking the Bayesian optimal design is that it needs heavy computation in the numerical calculations. To overcome this problem, a clustering method can be applied. Previous researches commonly used the k-means clustering method to alternate the full posterior distribution (Dror & Steinberg, 2006; McCallum & Bornkamp, 2015). K-means clustering reduces the computation of the high-dimensional data by partitioning $N$ observations into $k$ clusters by minimizing the within-cluster sum of squares (Hartigan & Wong, 1979). Various of other clustering methods have been developed such as kernel k-means and fuzzy c-means. Each of them has their own features. For instance, the kernel k-means handles the non-linear structure and the fuzzy c-means allows a data point to belong to two or more clusters (Dhillon et al., 2004; Welling, 2013; Dunn, 1973; Bezdek, 1981). In this research, we are also interested in comparing the performance of the three clustering methods to see if different clustering methods change the performance of the Bayesian c-optimal design. Without loss of generality, c-optimal design constructed in this paper is locally c-optimal design.

Section 2 shows the statistical model along with the information matrix. In Section 3, adaptive Bayesian c-optimal design for estimating multiple $\text{ED}_{p}s$ which account for the parameter uncertainty are proposed. Several scenarios of simulations are set up to check the performance of the proposed adaptive Bayesian c-optimal design compared with other traditional designs such as uniform design and D-optimal design. Additionally, we compare our design with the two-stage robust c-optimal design that was proposed by Zhang and Hyun (2016). All the design performance are evaluated by estimating multiple $\text{ED}_{p}s$ over a wide range of mis-specified model parameters. The final discussion is given in Section 4.

2. Statistical model

It is often observed that the dose-response relationship follows a sigmoid curve. It leads to the 4PL model which is frequently used for many toxicological and biological systems. It is a highly flexible nonlinear model taking into account the minimum dose, maximum dose, $\text{ED}_{50}$ , and the slope of the dose-response curve. The mean response $\eta$ for the 4PL model at a given dose ${x}_{i}$ is defined as

$\displaystyle\eta({x}_{i},\theta)={\theta}_{1}+({\theta}_{2}-{\theta}_{1})% \frac{{{X}_{i}}^{{\theta}_{4}}}{{{X}_{i}}^{{\theta}_{4}}+{{\theta}_{3}}^{{% \theta}_{4}}}$ (1)

where ${\theta}_{1}$ is the minimum mean effect; ${\theta}_{2}$ is the maximum mean effect; ${\theta}_{3}$ is the $\text{ED}_{50}$ , and ${\theta}_{4}$ is the slope parameter.

Let $Y$ be a continuous response at a given dose range $x=[0,{D}_{\max}]$ , where 0 is the placebo dose and ${D}_{\max}$ is the maximum dose. The continuous response ${Y}_{ij}$ for the $i^{\rm th}$ dose level that applies to $j^{\rm th}$ subject can be written as

$\displaystyle{Y}_{ij}=\eta({x}_{i},\theta)+{\varepsilon}_{ij},{\varepsilon}_{% ij}\sim N(0,{\sigma}^{2})$ (2)

where $j=1,2,3,\ldots,{n}_{i}$ , $i=1,2,3,\ldots,k$ ; ${\varepsilon}_{ij}$ is the residual error with the assumption of independent and normally distributed with variance ${\sigma}^{2}$ . The total sample size is $N={n}_{1}+\dots+{n}_{i}$ ; Let $\xi={\{({x}_{i},{w}_{i})\}}^{k}_{1}$ denotes the approximate design that allocates $N{w}_{i}$ subjects into ${x}_{i}$ dose level. Then the normalized Fisher information matrix under the 4PL model is

$\displaystyle M(\xi;\theta)=\frac{1}{\sigma^{2}}\sum^{k}_{i=1}{{w}_{i}f({x}_{i% },\theta){f({x}_{i},\theta)}^{T}}$ (3)

where $f({x}_{i},\theta)$ is the partial derivative of the mean response $\eta({x}_{i},\theta)$ with respect to the model parameters $\theta$ (matrix details see Padmanabhan & Dragalin, 2010; Zhang & Hyun, 2016). The Fisher information matrix $M(\xi;\theta)$ contains the design and model parameters. Hence, it is a key component in searching optimal designs in the next Section.

3. Adaptive bayesian c-optimal designs

3.1 Bayesian c-optimality

The explicit expression of the target dose $\text{ED}_{p}$ is ${\theta}_{3}$ ${(\frac{p}{1-p})}^{\frac{1}{{\theta}_{4}}}$ , where $p$ represents the $p$ % of the maximum treatment response. Once $p$ is fixed, the $\text{ED}_{p}$ becomes a function of model parameters, and it only depends on ${\theta}_{3}$ and ${\theta}_{4}$ . To estimate the function of model parameters, c-optimality can be applied, which minimizes the variance of estimating parameters. Let $\widehat{\text{ED}_{p}}$ denotes the maximum likelihood estimate of the $\text{ED}_{p}$ ; then the c-optimality can be redefined as minimizing the variance of estimating the $\text{ED}_{p}$ ,

$\displaystyle{\Psi}_{c}=\min_{\xi}\{Var(\widehat{\text{ED}_{p}})\}={\min_{\xi}% \{[\text{ED}_{p}^{\prime}]}^{T}M(\xi;\theta)^{-}\text{ED}_{p}^{\prime}\}$ (4)

where $\text{ED}_{p}^{\prime}$ is partial derivative of the $\text{ED}_{p}$ with respect to $\theta$ ; ${M(\xi;\theta)}^{-}$ is the generalized inverse of the Fisher information matrix in Eq. (3).

Our research interest is conducting one single study to estimate multiple $\text{ED}_{p}$ s taking into account the parameter uncertainty. One method to address the parameter uncertainty is applying the Bayesian paradigm at the second stage. Instead of using a single set of the parameters, the Bayesian approach considers the design incorporated with a prior probability for the parameters. To estimate the target dose $\text{ED}_{p}$ , we adopt the Bayesian c-optimality which minimizes the weighted average of the variance for estimating the $\text{ED}_{p}$ (Atkinson & Donev, 1992):

$\displaystyle{\Psi}_{Ba}=\min_{\xi}\left\{\sum^{m}_{i=1}P\left({\theta}_{i}% \right)Var\left(\widehat{\text{ED}_{p}}|{\theta}_{i}\right)\right\}$ (5)

where $P({\theta}_{i})$ denotes the prior probability given to ${\theta}_{i}$ , and $\sum^{m}_{i=1}{P({\theta}_{i})}=1$ . If there is no prior knowledge of the probability $P({\theta}_{i})$ and all values of the parameters are considered equally important, we can simply set $P\left({\theta}_{i}\right)=\frac{1}{m}$ . The fact is that the Bayesian optimal design truly depends on the prior distribution of the parameters. To reduce the prior dependency, the posterior distribution $P({\theta}_{i}|y)$ is applied instead, which can be obtained bases on the response data collected in the first stage design and the pre-determined prior density (McCallum & Bornkamp, 2015).

Two stages are involved in searching the adaptive Bayesian c-optimal design: In the first stage, a small proportion of the subjects is assigned to a fixed design such as uniform design. In the second stage, the posterior distribution of the model parameters is obtained based on the prior distribution and the information that is collected from the first stage. Next, the generated posterior distribution is incorporated into the Bayesian c-optimality.

There are two challenges to construct adaptive Bayesian c-optimal design at the second stage: (1) computing the undefined posterior distribution; (2) heavy numerical evaluations of the optimal criterion. In order to consider conservative performance of the adaptive Bayesian c-optimal design, the Jefferys prior distribution is used: $P(\theta)={\left|M({\xi}_{1};\theta)\right|}^{1/2}$ , where ${\xi}_{1}$ is the fixed design at the first stage. Due to the Jeffreys prior distribution is included, the posteriors distribution becomes an undefined distribution. In this paper, we simulate sampling draws by using the Markov Chain Monte Carlo (MCMC) with Metropolis-Hastings (M-H) algorithm to generate 10,000 samples with 1000 burn-in from the posterior distribution of the parameters. As an illustration purpose, we adopt the same experimental setup values from Padmanabhan and Dragalin (2010) to perform the numerical search in the first stage design. That is, the dose range $x=[0,8]$ and the initial parameters $\theta=(0,-1.7,4,5)$ . Under the 4PL model, let ${\sigma}^{2}=\sqrt{2}$ , and we simply assign a small sample size of 40 subjects into a fixed design, say four points uniform design,

$\displaystyle\xi_{1}=\begin{pmatrix}0.001,&2.67,&5.33,&8\\ 0.25,&0.25,&0.25,&0.25\end{pmatrix},$

within the continuous dose range from 0 to 8. Once all the initial values are fixed, we can start to generate the response data and use for searching the posterior distribution of the parameters. Since the c-optimal design under 4PL does not depend on ${\theta}_{1}$ and ${\theta}_{2}$ , we focus on generating the posterior distributions for ${\theta}_{3}$ and ${\theta}_{4}$ . Based on the collected information from the first stage, a sample size of 9000 for ${\theta}_{3}$ and ${\theta}_{4}$ are generated, respectively. Their histogram and density plots are given in Fig. 1. The adaptive Bayesian c-optimal design using the sampling parameters from the posterior distribution are obtained in the following subsection.

Figure 1.

Histogram and density plots of the posterior distributions for ${\theta}_{{3}}$ and ${\theta}_{{4}}{}$ with generated sample size of 10,000 and 1000 burn-in, respectively.

3.2 Adaptive Bayesian c-optimal design for estimating the

{\text{ED}}_{{50}}

For simplicity, we start with searching for an adaptive Bayesian c-optimal design for estimating the $\text{ED}_{50}$ , which minimizes the weighted average of the variance for estimating the $\text{ED}_{50}$ . The directional derivative of the Bayesian c-optimality criterion in Eq. (5) leads to the sensitivity function as below:

$\displaystyle\sum^{m}_{i=1}{{\lambda}_{i}\frac{{\left\{{f}^{T}\left(x,{\theta}% _{i}\right){M\left({{\xi}_{Ba}};{\theta}_{i}\right)}^{-}\text{ED}^{\prime}_{50% }\right\}}^{2}}{{{\left[\text{ED}^{\prime}_{50}\right]}^{T}M\left({{\xi}_{Ba}}% ;{\theta}_{i}\right)}^{-}\text{ED}^{\prime}_{50}}}\leqslant 1$ (6)

where ${\lambda}_{i}$ is the weight of each set of the ${\theta}_{i}$ and is set to $\frac{1}{m}$ with $m$ sampling draws. According to the General Equivalence Theorem, the equal sign holds, if and only if, $x$ is one of the design points in the adaptive Bayesian c-optimal design ${\xi}_{Ba}$ . One issue for the design implementation is that 9000 parameter samples from the posterior distribution requires a computation of 9000 iterations of the Fisher information matrix for searching one design. It is very time-consuming to complete the whole procedure.

In illustration of using the full posterior distribution, we randomly select a sample size of 200 from the 9000 samples to obtain the adaptive Bayesian c-optimal design for estimating the $\text{ED}_{50}$ . The Bayesian optimal design is searched by adopting the modified Yang-Biedermann-Tang (YBT) algorithm which proposed by Hyun and Wong (2015). Instead of starting with some arbitrary design as an initial design, the modified YBT algorithm starts with a better initial design that suggested by running a traditional algorithm for several iterations (such as V-algorithm). At each iteration, the modified YBT algorithm adds the design point that maximizes the sensitivity function into the current design to create a new design. The weights for the new design point are computed using the Newton-Raphson method by optimizing the design criterion over the design space (Yang et al., 2013; Hyun & Wong, 2015; Zhang & Hyun, 2016). The obtained Bayesian c-optimal design is verified by the Equivalence Theorem in Fig. 2. The adaptive Bayesian c-optimal design, ${\xi}_{Ba}$ , for estimating the $\text{ED}_{50}$ by using the 200 sample parameters from the posterior distribution is given by:

$\displaystyle{\xi}_{Ba}=\begin{pmatrix}0.0001,&0.01,&0.42,&1.28,&1.68,&2.20,&3% .86,&7.99\\ 0.02,&0.03,&0.12,&0.03,&0.07,&0.20,&0.31,&0.22\end{pmatrix}$ (7)

One can observe that the adaptive Bayesian c-optimal design for estimating the $\text{ED}_{50}$ contains 8 design points including the lower bound and the upper bound of dose levels. Dose levels 0.42, 2.20, 3.86, and 7.99 take most of the design weights than others. It suggests that most of the subjects are allocated into these four dose levels.

Figure 2.

Plot of the sensitivity function of adaptive Bayesian c-optimal design for estimating the $\text{ED}_{50}$ under the full posterior distribution with 200 sets of parameters that randomly selected from the generated 9000 samples.

McCallum and Bornkamp (2015) proposed an efficient alternative to reduce the laborious evaluations by adopting a k-means clustering method instead of using full posterior samples. They also proved that the k-means clustering provided good theoretically approximations when number of clusters $k=10$ . In our study, we classify $m=9000$ samples into $k=10$ clusters using k-means clustering algorithm to conduct the adaptive Bayesian c-optimal design for estimating the $\text{ED}_{50}$ . In this scenario, the ${\lambda}_{i}$ in the sensitivity function (6) becomes the weight of each clusters, which is the proportion of the samples associate with each cluster center. As a comparison, kernel k-means and fuzzy c-means clustering methods are studied in searching the Bayesian design. Constantly, the kernel k-means and fuzzy c-means also use $k=10$ clusters. K-means, kernel k-means, and fuzzy c-means clustering can be easily computed in R software.

The adaptive Bayesian c-optimal design for estimating the $\text{ED}_{50}$ using k-means, kernel k-means, and fuzzy c-means clustering methods are given in Table 1. All the searched designs are verified by the General Equivalence Theorem (verification plots are not shown here). Table 1 states that the design points and design weights are varied from the three different Bayesian designs, but all the designs contain 8 design points include lower bound and upper bound of the dose range.

For further investigation of the obtained designs performance, we compute the design efficiencies for the cases using the full posterior, k-means, kernel k-means, and fuzzy c-means. Design efficiency measures how a design performs with respect to some optimality criterion. Our research goal is to estimate the target dose $\text{ED}_{p}$ , so we use c-efficiency which is the ratio of the variance of a design $\xi$ against the variance of c-optimal design ${{\xi}^{*}}_{{ED}_{p}}$ for estimating the $\text{ED}_{p}$ . It is defined a

$\displaystyle\textit{Eff}_{\text{ED}_{p}}(\xi)=\frac{[\text{ED}^{\prime}_{p}]^% {T}M({\xi}^{*}_{\text{ED}_{p}};\Theta)^{-}\text{ED}^{\prime}_{p}}{[\text{ED}^{% \prime}_{p}]^{T}M(\xi;\Theta)^{-}\text{ED}^{\prime}_{p}}$ (8)

Table 1

Adaptive Bayesian c-optimal designs that obtained by using k-means, kernel k-means, and fuzzy c-means clustering methods. ${\xi}_{Ba.k}$ , ${\xi}_{Ba.kk}$ , ${\xi}_{Ba.fc}$ represents the adaptive Bayesian c-optimal design with k-means clustering, kernel k-means clustering, and fuzzy c-means clustering, respectively.

Clustering methods	$\xi$
K-means	${\xi}_{Ba.k}=\begin{pmatrix}0.0001,&0.98,&1.46,&2.27,&2.49,&4.05,&5.40,&7.99\\ 0.11,&0.06,&0.03,&0.20,&0.05,&0.24,&0.10,&0.22\\ \end{pmatrix}$
Kernel k-means	${\xi}_{Ba.kk}=\begin{pmatrix}0.0001,&0.04,&0.97,&1.64,&2.16,&2.61,&4.49,&7.99% \\ 0.09,&0.04,&0.02,&0.09,&0.07,&0.18,&0.28,&0.22\\ \end{pmatrix}$
Fuzzy c-means	${\xi}_{Ba.fc}=\begin{pmatrix}0.0001,&0.54,&1.49,&2.46,&2.67,&4.34,&4.63,&7.99% \\ 0.11,&0.02,&0.12,&0.12,&0.10,&0.22,&0.07,&0.23\\ \end{pmatrix}$

The locally c-optimal design ${{\xi}^{*}}_{\text{ED}_{p}}$ for estimating the $\text{ED}_{p}$ provides the minimum variance, thus the c-efficiency is always between 0 and 1. Higher value of c-efficiency indicates better performance of the design $\xi$ .

As a comparison, traditional designs, such as uniform design and D-optimal design, are studied. D-optimal design is widely used for estimating the dose response curve. Even though it is not the goal of D-optimal design, it is also commonly used as an illustration for estimating a particular dose such as $\text{ED}_{50}$ or $\text{ED}_{90}$ (Dragalin et al., 2007; Li et al., 2008). The searched D-optimal design for estimating $\text{ED}_{50}$ using D-optimality is given by

$\displaystyle{\xi}_{D}=\begin{pmatrix}0.001,&3.17,&4.71,&7.92\\ 0.25,&0.25,&0.25,&0.25\end{pmatrix}.$

Considering the obtained adaptive Bayesian c-optimal designs all contain 8 design points, we compare with the 8 points uniform design,

$\displaystyle{\xi}_{U1}=\begin{pmatrix}0.001,&1.14,&2.29,&3.43,&4.57,&5.71,&6.% 86,&8.00\\ \frac{1}{8},&\frac{1}{8},&\frac{1}{8},&\frac{1}{8},&\frac{1}{8},&\frac{1}{8},&% \frac{1}{8},&\frac{1}{8};\end{pmatrix}$

Additionally, we are interested in examining how the Bayesian technique enhances the performance of traditional uniform design by searching Bayesian optimal weights for each uniform design points. To find the Bayesian optimal weights, we maximize the Bayesian c-optimal criterion function over the weights for the given 8 design points. The searched 8 points uniform design with Bayesian optimal weights for estimating $\text{ED}_{50}$ is

$\displaystyle{\xi}_{U2}=\begin{pmatrix}0.001,&1.14,&2.28,&3.43,&4.57,&5.71,&6.% 86,&8.00\\ 0.11,&0.09,&0.23,&0.08,&0.24,&0.01,&0.00,&0.24\end{pmatrix}.$

Note that after applying Bayesian technique, the design weights are not equally distributed into the 8 design points in the traditional way. The Bayesian optimal weights suggest to allocate most of subjects into dose levels into lower bound 0.001, middle levels 2.28 and 4.57, and upper bound 8.00. We are also interested in checking how the two-stage robust c-optimal design performs under the wide range of 9000 sampling parameters, which was proposed by Zhang and Hyun (2016). The summary statistics of 25% quartile, median, mean, and 75% quartile of the c-efficiencies are provided in Table 2.

From the summary statistics table, we note that there are very slight differences among the four adaptive Bayesian c-optimal designs and the uniform design with Bayesian optimal weights when estimate the $\text{ED}_{50}$ . One potential reason that the uniform design with Bayesian optimal weights performs as well as the Bayesian optimal designs could be that it contains 8 equally spaced design points which is very close to the Bayesian optimal designs. If we reduce the 8 uniform design points to 4 or 5 design points, the uniform design with Bayesian optimal weights would perform much poorer.

Table 2

Summary statistics of c-efficiencies of the designs for estimating the $\text{ED}_{50}$ under the 9000 $\theta s$ . ${\xi}_{Ba}$ , ${\xi}_{Ba.k}$ , ${\xi}_{Ba.kk}$ , and ${\xi}_{Ba.fc}$ represent the adaptive Bayesian c-optimal design with full posterior distribution, k-means clustering, kernel k-means clustering, and fuzzy c-means clustering, respectively; ${\xi}_{Two}$ represents the two-stage c-optimal robust design; ${\xi}_{D}$ represents the D-optimal design. ${\xi}_{U1}$ represents the 8 points uniform design; ${\xi}_{U2}$ represents the 8 points uniform design with Bayesian optimal weights

$\xi$	25% Q	Median	Mean	75% Q
${\xi}_{Ba}$	0.4981	0.5931	0.5763	0.6601
${\xi}_{Ba.k}$	0.4996	0.5763	0.5596	0.6189
${\xi}_{Ba.kk}$	0.5074	0.5704	0.5696	0.6327
${\xi}_{Ba.fc}$	0.5062	0.5777	0.5723	0.6438
${\xi}_{Two}$	0.3219	0.4762	0.5022	0.6662
${\xi}_{D}$	0.1388	0.3855	0.3649	0.5890
${\xi}_{U1}$	0.4481	0.4775	0.4739	0.5026
${\xi}_{U2}$	0.5018	0.5800	0.5668	0.6382

As expected, the two-stage robust c-optimal design for estimating the $\text{ED}_{50}$ works worse compared with the Bayesian optimal designs. The mean efficiency of the two-stage robust c-optimal design is around 7% lower and its lower quartile is much lower. This implies the two-stage robust c-optimal design is not as robust as the adaptive Bayesian c-optimal design for mis-specified parameter values. Table 2 shows that the median and the mean design efficiencies for the Bayesian designs and uniform design with the Bayesian optimal weights are around 10% higher than the traditional uniform design, and 20% higher than the D-optimal design. Overall, we conclude that after reducing the parameter dependency, the adaptive Bayesian c-optimal designs still work better than the uniform design, D-optimal design and the two-stage robust c-optimal design. Since D-optimal design works the worst under a wide range of mis-specified parameter values, we exclude D-optimal design in our following design comparisons.

3.3 Adaptive Bayesian c-optimal design for estimating multiple

{\text{ED}}_{{p}}

From the previous section, we show that k-means, kernel k-means, and fuzzy c-means perform similarly when compared with the full posterior method. For simplicity, we use the k-means algorithm to extend our study in searching the adaptive Bayesian c-optimal design for estimating multiple $\text{ED}_{p}$ s. The sensitivity function can be updated from the Eq. (6) by summarizing the weighted multiple values of $\text{ED}_{p}$ :

$\displaystyle\sum^{k}_{i=1}{{\lambda}_{i}\sum^{l}_{j=1}{{w}_{j}\frac{\{{f}^{T}% (x,{\theta}_{i}){M({{\xi}_{Ba}};{\theta}_{i})}^{-}\text{ED}^{\prime}_{j}\}^{2}% }{{[\text{ED}^{\prime}_{j}]}^{T}M({{\xi}_{Ba}};{\theta}_{i})^{-}\text{ED}^{% \prime}_{j}}}}\leqslant 1$ (9)

where $l$ is the number of $\text{ED}_{p}$ s, $k$ is the number of the clusters, ${\lambda}_{i}$ represents the weight of each cluster, and ${w}_{j}$ represents the weight of each $\text{ED}_{p}$ . Here we set $l=5$ and $k=10$ . Note that the values of $l,k,\lambda$ , and $w$ in the Eq. (9) can be specified as any value depending on the research interests. Again, by the Equivalence Theorem, the equal sign holds if, and only if, $x$ is one of the dose levels in the adaptive Bayesian c-optimal design, ${\xi}_{Ba}$ . Again, the modified YBT algorithm is applied to search the design for estimating multiple $\text{ED}_{p}$ s $(\text{ED}_{10},\text{ED}_{30},\text{ED}_{50},\text{ED}_{70},\text{ED}_{90})$ , ${\xi}_{Ba.\text{ED}_{p}s}$ as below:

$\displaystyle{\xi}_{Ba.\text{ED}_{p}s}=\begin{pmatrix}0.0001,&0.23,&1.04,&1.6,% &2.56,&2.75,&3.7,&5.03,&5.25,&7.99\\ 0.08,&0.07,&0.04,&0.13,&0.16,&0.03,&0.15,&0.13,&0.03,&0.18\end{pmatrix}$ (10)

Again, the ${\xi}_{Ba.\text{ED}_{p}s}$ is verified by General Equivalence Theorem.

Note that we illustrate our study by including the five $\text{ED}_{p}$ s above to conduct the adaptive Bayesian c-optimal design, but it can be any $\text{ED}_{p}$ s values. Next, we investigate the performance of the design for estimating other $\text{ED}_{p}$ s. Considering cover all the possible values of $\text{ED}_{p},$ we use ten $\text{ED}_{p}$ s ( $\text{ED}_{10},\text{ED}_{20},\text{ED}_{30},\text{ED}_{40},$ $\text{ED}_{50}$ , $\text{ED}_{60}$ , $\text{ED}_{70}$ , $\text{ED}_{80}$ , $\text{ED}_{90}$ , and $\text{ED}_{99}$ ) to compute the design efficiencies. We randomly select ten sets of model parameters from the 9000 samples and check how the adaptive Bayesian c-optimal design works to estimate the ten $\text{ED}_{p}s$ under the ten sets of parameter values. The randomly selected ten sets of model parameters are shown in Table 3. In the c-efficiency computation, the adaptive Bayesian c-optimal design, the uniform designs, and the two-stage robust c-optimal design are all fixed design. Take estimating the $\text{ED}_{10}$ under the model parameter ${\theta}_{(1)}=(0,-1.7,2.6,3.6)$ as an example, first we search the c-optimal design for estimating the $\text{ED}_{10}$ under ${\theta}_{(1)}$ , then we compute the c-efficiency by using one of the fixed designs, for instance, the adaptive Bayesian c-optimal design ${\xi}_{Ba.\text{ED}_{p}s}$ , against the c-optimal design. For each design ( ${\xi}_{Ba.\text{ED}_{p}s}$ , ${\xi}_{Two}$ , ${\xi}_{U1}$ , ${\xi}_{U2}$ ), this procedure repeats 9000 times until all the c-efficiencies are computed based on the 9000 parameter samples. The detailed summary statistics of the c-efficacies are shown in Table 4.

Table 3

Ten sets of model parameters that are randomly selected from the generated posterior distribution

$\theta s$
${\theta}_{(1)}=(0,1.7,2.6,3.6)$	${\theta}_{(6)}=(0,1.7,2.4,2.8)$
${\theta}_{(2)}=(0,1.7,2.9,4.8)$	${\theta}_{(7)}=(0,1.7,2.8,4.4)$
${\theta}_{(3)}=(0,1.7,2.5,4.0)$	${\theta}_{(8)}=(0,1.7,2.9,4.5)$
${\theta}_{(4)}=(0,1.7,3.5,3.1)$	${\theta}_{(9)}=(0,1.7,2.1,4.6)$
${\theta}_{(5)}=(0,1.7,2.3,4.6)$	${\theta}_{(10)}=(0,1.7,3.63.2)$

Table 4 shows that when estimate the $\text{ED}_{10}$ and the $\text{ED}_{30}$ , the mean efficiencies of the two-stage robust c-optimal design are slightly higher than the other designs in the simulation. However, for estimating the $\text{ED}_{50}$ , $\text{ED}_{70}$ , and $\text{ED}_{90}$ , adaptive Bayesian c-optimal design and the 8 points uniform design with Bayesian optimal weights work better than the two-stage robust c-optimal design. This suggests that two-stage robust c-optimal design sometimes works better, but sometimes not when compared with the adaptive Bayesian c-optimal design. No substantial difference is observed between the adaptive Bayesian c-optimal design and 8 points uniform design with Bayesian optimal weights, which further support our discussion in the previous sections. However, when we further investigate the two designs, the 25% quartile of the adaptive Bayesian c-optimal design are always higher than the Bayesian enhanced 8 points uniform design. It suggests that the adaptive Bayesian c-optimal design is much safer to use when estimates a wide range of the $\text{ED}_{p}s$ since its worst case is always better than others. As we expected, the traditional 8 points uniform design always works the worst with varied values of the $\text{ED}_{p}$ and the mis-specified model parameters.

Table 4

Summary statistics of c-efficiencies of the designs for estimating multiple $\text{ED}_{p}$ s ( $\text{ED}_{10},\text{ED}_{30},\text{ED}_{50},$ $\text{ED}_{70},$ and $\text{ED}_{90}$ ). ${\xi}_{Ba.\text{ED}_{p}s}$ represents the adaptive Bayesian c-optimal design; ${\xi}_{Two}$ represents the two-stage robust c-optimal design; ${\xi}_{U1}$ represents the traditional 8 points uniform design; ${\xi}_{U2}$ represents the 8 points uniform design with Bayesian optimal weights

$\text{ED}_{p}$	$\xi$	25% Q	Median	Mean	75% Q
$\text{ED}_{10}$	${\xi}_{Ba.\text{ED}_{p}s}$	0.4030	0.4657	0.4656	0.5260
	${\xi}_{Two}$	0.3536	0.4835	0.4837	0.5817
	${\xi}_{U1}$	0.3303	0.3763	0.3884	0.4654
	${\xi}_{U2}$	0.3656	0.4387	0.4391	0.5111
$\text{ED}_{30}$	${\xi}_{Ba.\text{ED}_{p}s}$	0.4715	0.5369	0.5203	0.5879
	${\xi}_{Two}$	0.4330	0.5579	0.5444	0.6621
	${\xi}_{U1}$	0.3855	0.4671	0.4411	0.5122
	${\xi}_{U2}$	0.4381	0.5349	0.5114	0.6049
$\text{ED}_{50}$	${\xi}_{Ba.\text{ED}_{p}s}$	0.5216	0.5550	0.5526	0.5830
	${\xi}_{Two}$	0.4316	0.5540	0.5268	0.6511
	${\xi}_{U1}$	0.4489	0.4783	0.4777	0.5032
	${\xi}_{U2}$	0.5138	0.5666	0.5561	0.6074
$\text{ED}_{70}$	${\xi}_{Ba.\text{ED}_{p}s}$	0.4916	0.5269	0.5314	0.5684
	${\xi}_{Two}$	0.3986	0.4907	0.4707	0.5396
	${\xi}_{U1}$	0.4177	0.4534	0.4596	0.5049
	${\xi}_{U2}$	0.4824	0.5375	0.5302	0.5820
$\text{ED}_{90}$	${\xi}_{Ba}$	0.5059	0.5570	0.5553	0.6155
	${\xi}_{Two}$	0.2969	0.4342	0.4326	0.5489
	${\xi}_{U1}$	0.4350	0.4918	0.4766	0.5469
	${\xi}_{U2}$	0.4748	0.5737	0.5418	0.6204

4. Discussion

A key objective of dose-finding trials is often to study the dose-response curve or to estimate target dose levels of interest such as $\text{ED}_{p}$ . In this research, our interest lies in conducting an optimal design to estimate multiple target doses, $\text{ED}_{p}s$ , simultaneously, whilst accounting for parameter dependency. We employ the two-stage design procedure to construct the adaptive Bayesian c-optimal design that incorporated the Bayesian paradigm in the second stage. It is shown that the proposed design reduces the parameter dependency without lose much design efficiency.

One challenge using the full posterior distribution to carry out the Bayesian paradigm in the second stage is the laborious evaluations of the information matrix for searching one Bayesian design. To overcome the heavy computation issue, three clustering methods such as k-means, kernel k-means, and fuzzy c-means are utilized as alternatives in constructing the Bayesian optimal design. Summary statistics of the simulation demonstrates that there is not much change among the performance of the three clustering methods Bayesian optimal designs.

From the c-efficiency histogram plots over the varied sets of parameter values, we state that the adaptive Bayesian c-optimal design and the uniform design with Bayesian optimal weights work reasonablely well for estimating varied values of $\text{ED}_{p}$ s compared with traditional uniform design and two-stage robust c-optimal design. However, a minimal difference exists between the adaptive Bayesian c-optimal design and the 8 points uniform design with Bayesian optimal weights. One potential reason as discussed earlier could be the uniform design with optimal weights has 8 equal spaced design points which are relatively close to the adaptive Bayesian c-optimal design. Comparable results are found in the later computation of the c-efficiencies. The c-efficiency results demonstrate that the two-stage robust c-optimal design works slightly better when estimating lower value of $\text{ED}_{p}$ ( $\text{ED}_{10}$ and $\text{ED}_{30}$ ); However, when estimating higher value of $\text{ED}_{p}$ ( $\text{ED}_{50},\text{ED}_{70}$ , and $\text{ED}_{90}$ ), adaptive Bayesian c-optimal design become better.

In summary, this research shows that the proposed adaptive Bayesian c-optimal design accounting for parameter uncertainty works fairly well for estimating multiple $\text{ED}_{p}$ s. Compared with the adaptive Bayesian c-optimal design, the two-stage robust c-optimal design sometimes works better, sometimes not. The adaptive Bayesian c-optimal design works similarly with the 8 points uniform design with Bayesian optimal weights. However, when estimates wide range of the $\text{ED}_{p}$ s, the Bayesian c-optimal design is much safer to use because its worst case is always better than the other designs. Another result of note is that the uniform design with optimal weights employing the Bayesian technique significantly enhanced the design performance compared with the traditional uniform design.

References

Albert

(2009). Bayesian Computation with R. Springer, New York.

Atkinson

A. C.

, & Donev

A. N.

(1992). Optimum Experimental Designs. Oxford University Press, Oxford.

Beloeil

Eurin

Thevenin

Benhamou

, & Mazoit

J. X.

(2007). Effective dose of nefopam in 80% of patients (ED80): A study using the continual reassessment method. British Journal of Clinical Pharmacology, 64, 686-693.

Bezdek

J. C.

(1981). Pattern Recognition with Fuzzy Objective Function Algoritms. Plenum Press, New York.

Chernoff

(1953). Locally optimal designs for estimating parameters. Annals of Mathematical Statistics, 39, 1266-1281.

Dette

(1996). A note on Bayesian c- and d-optimal design in nonlinear regression models. The Annals of Statistics, 24(3), 1225-1234.

Dhillon

I. S.

Guan

Y. Q.

, & Kulis

(2004). Kernel k-means, Spectral Clustering and Normalized Cut. Retrieved from http://www.cs.utexas.edu/users/inderjit/public_papers/kdd_spectral_kernelkmeans.pdf.

Dragalin

Hsuan

, & Padmanabhan

S. F.

(2007). Adaptive Designs for Dose-finding Studies Based on Sigmoid Emax Model. Journal of Biopharmaceutical Statistics, 17, 1051-1070.

Dror

, & Steinberg

D. M.

(2006). Robust experimental design for multivariate generalized linear models. Technometrics, 48(4), 520-529.

10.

Dunn

J. C.

(1973). A fuzzy relative of the ISODATA process and its use in detecting compact well-separated clusters. Journal of Cybernetics, 3, 32-57.

11.

Hartigan

J. A.

, & Wong

M. A.

(1979). A k-means clustering algorithm. Applied Statistics, 28, 100-108.

12.

Hyun

S. W.

, & Wong

W. K.

(2015). Multiple-objective optimal designs for studying the dose response function and interesting dose levels. The International Journal of Biostatistics, 11(2), 1557-4679.

13.

Louie

Kaw

Liu

Jumbe

Miller

M. H.

, & Drusano

G. L.

(2001). Pharmacodynamics of daptomycin in a murine thigh model of staphylococcus aureus infection. Antimicrobial, Agents and Chemotherapy, 45(3), 845-851.

14.

, & Majumdar

(2008). D-optimal design for logistic models with three and four parameters. Journal of Statistical Planning and Inference, 138, 1950-1959.

15.

McCallum

, & Bornkamp

(2015). Accounting for parameter uncertainty in two-stage designs for Phase II dose-response studies. Modern Adaptive Randomized Clinical Trials: Statistical and Practical Aspects, Boca Raton: Chapman & Hall/CRC Press, 427-450.

16.

Miller

Guilbaud

, & Dette

(2007). Optimal designs for estimating the interesting part of a dose-effect curve. Journal of Biopharmaceutical Statistics, 17, 1097-1115.

17.

Padmanabhan

S. K.

, & Dragalin

(2010). Adaptive dc-optimal design for dose finding based on a continuous efficacy endpoint. Biometrical Journal, 52, 836-52.

18.

Ting

(2006). Dose Finding in Drug Development. Springer, New York.

19.

Wang

T. H.

, & Yang

(2014). Adaptive optimal designs for dose-finding studies based on sigmoid Emax models. Journal of Statistical Planning and Inference, 122, 188-197.

20.

Welling

(2013). Kernel k-means and spectral clustering. Retrieved from http://www.ics.uci.edu/∼welling/teaching/273ASpring09/.SpectralClustering.pdf.

21.

Wong

W. K.

(1994). Comparing robust properties of a, d, e and g-optimal designs. Computational Statistics & Data Analysis, 18, 441-448.

22.

Yang

Biedermann

, & Tang

(2013). On optimal deigns for nonlinear models: A general and efficient algorithm. Journal of the American Statistical Association, 108, 1411-1420.

23.

Zhang

, & Hyun

S. W.

(2016). Robust c-optimal design for estimating multiple e⁢dps under the 4-parameter logistic model. Statistics, Optimization and Information Commutating, 4, 278-288.