Spline-based procedures for dose-finding studies with active control.

In a dose-finding study with an active control, several doses of a new drug are compared with an established drug (the so-called active control). One goal of such studies is to characterize the dose-response relationship and to find the smallest target dose concentration d(*), which leads to the same efficacy as the active control. For this purpose, the intersection point of the mean dose-response function with the expected efficacy of the active control has to be estimated. The focus of this paper is a cubic spline-based method for deriving an estimator of the target dose without assuming a specific dose-response function. Furthermore, the construction of a spline-based bootstrap CI is described. Estimator and CI are compared with other flexible and parametric methods such as linear spline interpolation as well as maximum likelihood regression in simulation studies motivated by a real clinical trial. Also, design considerations for the cubic spline approach with focus on bias minimization are presented. Although the spline-based point estimator can be biased, designs can be chosen to minimize and reasonably limit the maximum absolute bias. Furthermore, the coverage probability of the cubic spline approach is satisfactory, especially for bias minimal designs.

1 Introduction

The aim of phase II dose-finding studies is the characterization of the dose–response relationship of the investigated drug. Sometimes, it is of interest to estimate a target dose, which leads to a desired response. An overview over the variety of approaches on how to deal with the planning and analysis of this type of studies is given by Ruberg 1,2 and by Bornkamp et al. 3. More specifically, Bretz et al. 4 and Pinheiro et al. 5 proposed a methodology that combines formal hypothesis testing for dose–response with flexible modeling of the dose–response relationship. For a more detailed discussion of their MCPMod approach, we refer to 6,7. Estimating the smallest dose with a clinically relevant effect in a parametric model can be seen as a calibration problem or inverse regression problem. This requires estimating the value of an independent variable that yields an expected outcome of the dependent variable equal to a prespecified value. There is extensive literature on this kind of calibration problems, in quality control and dose estimation 8–11. Whereas these approaches focus on a target dose that is defined by a specific fixed effect difference to placebo, there is an increasing interest in describing the dose–response relative to an active comparator to be included in the dose-finding trial. This, however, has found little attention in the statistical literature with the notable exceptions of 12,13. In addition, the use of a placebo may be unethical in some situations. If the dose-finding study includes an active control, the target dose d* is without loss of generality defined as the smallest dose that leads to the same expected efficacy as the active comparator. In this paper, the focus will lie on phase II dose-finding studies with active control where the underlying dose–response function is unknown, and no assumptions regarding the shape are made in the analysis as well as estimating the target dose d* and the construction of a corresponding CI.

We found little in statistical literature addressing this kind of problem except the approach of Kirby et al. 14, which uses a cubic smoothing spline, with cross-validation for the smoothing parameter to estimate the dose–response function but does not include an active comparator in the statistical model. Additionally, the construction of a CI for the target dose was not investigated. Further, the approach of Dilleen et al. 15, which includes an active control, only assumes a monotone response function and uses a linear spline interpolation to estimate the target dose and therefore fits the problem more precisely. In order to estimate the target dose without assumptions on the underlying dose–response shape, a cubic spline approach 16 will be used, and a cubic spline-based bootstrap method for constructing a CI 17 for the target dose will be presented. This approach will be compared by simulation studies with the linear spline approach from Dilleen et al. 15 as well as two parametric approaches, which use an approximation of the standard error by the Δ-method to construct corresponding CIs 18. The spline-based point estimators of the target dose are characterized in terms of bias, whereas coverage probabilities are presented for the 95% CI for all methods. Additionally, the possible benefits of a smoothed cubic spline over the cubic spline interpolation are investigated in a small simulation study. Further, design considerations for the introduced spline approaches with the focus on bias minimization will be presented in the context of dose-finding studies with active control.

Typically, phase II dose-finding studies include three to five doses of an experimental drug, a dose of an active comparator, and possibly a placebo control. The range of the used total sample sizes varies from 100 to more than 1000 subjects per study. Examples include the study of Krum et al. 19 investigating the effect of an endothelin receptor antagonist called bosentan on blood pressure in patients with essential hypertension. In this study, 293 patients were randomly assigned to receive placebo or one of the four oral doses of bosentan (100, 500, 1000, and 2000 mg per day) or the angiotensin-converting enzyme inhibitor enalapril as active control for 4 weeks, as well as a study presented by Chapple et al. 20, which investigated the effect of solifenacin on patients with symptomatic idiopathic detrusor overactivity. A total of 255 patients were randomized to receive placebo or one of the four doses of solifenacin (2.5, 5, 10, or 20 mg per day) or tolterodine as active control. The primary endpoints were voids/24 h and mean volume voided. The motivating example used in this paper is a study by Nauck et al. 21 in type 2 diabetes. The study was a double-blinded, double-dummy, placebo-controlled and active-controlled, parallel-group design over a period of 26 weeks. A total number of N = 1091 subjects, which were randomly assigned (2:2:2:1:2) to liraglutide (0.6, 1.2, or 1.8 mg/day) once daily, to placebo or to the active control glimepiride were enrolled. All treatments were in combination therapy with metformin (1 g twice daily). The primary endpoint was defined by the change in the percentage of glycated hemoglobin in the blood HbA1c at the end of the study. The secondary endpoints included the change in fasting plasma glucose (FPG) and body weight from baseline to the end of the study. For the purpose of illustration, the decrease in FPG in millimole per liter (mmol/l) will be used as endpoint for our motivating example. In the study, the standard deviations (SD) of the FPG values are only presented at the beginning and at the end of the study. For the decrease in FPG, no SD was available from the publication, and therefore, we set the SD to 1.8 mmol/l for all dose levels as well for the active control. The study results with assumed standard error are displayed in Figure 1. These were used to motivate the simulation scenarios.

Figure 1

FPG decrease of the experimental drug dose levels and active control at the end of the study displayed as mean± assumed standard error (SE) reported in 21.

This paper is organized as follows. In section 1, notation and statistical model are introduced. In section 2 and section 3, the properties of the cubic spline approach are investigated. More precisely, the estimation bias, the coverage probability of the CI as well as the median CI length are assessed in extensive simulations. Furthermore, the potential benefit of a smoothed spline is investigated. In Section 4, methods of deriving a bias-based optimal design for the use of cubic spline interpolation are presented. Finally, we close with a brief discussion and some recommendations in Section 5. All technical proofs can be found in Appendix A.

2 Spline-based estimation of the target dose

2.1 Notation and statistical model

In the following, bold letters such as or define a vector, and bold capital letters such as or a matrix, if not defined otherwise. Further denotes ′ the transpose of and −1 the inverse matrix. Especially, defines the identity matrix of dimension n, and 1 is the vector containing n times the element 1. The symbol ∼ defines ‘distributed as’ and defines ‘should be equal’.

For k groups and the dose levels d1person in the i-th dose level can be written as Y=f(d) + ε, with dose–response function f(d) and normally distributed error terms for i = 1,…,k,j = 1,…,n. The random variables of the active control (ac) are given by Y=μ + ε, with expected value μ and error terms for person j = 1,…,n. Let be the sample size of all dose levels and N = n+n the total sample size. Further describes the random vector of the k dose levels, the vector of the active control, and the vector of all random variables. The response function f(d) as well as the first, second and fourth derivatives and should exist and be continuous on the investigated dose range. No other assumptions will be made regarding the shape of f(d). To determine the target dose d*, the smallest dose that solves the equation f(d*) = μ must be found.

In general, the dose–response curve as well as the expected value of the active control μ are unknown and must be estimated from the data. In the next section, a method to construct a robust estimating procedure for the target dose without any additional assumptions on the dose–response function will be presented. Therefore, the term robust will be used if an estimator is robust against the different possible shapes of the function f(d) and does not depend on additional assumptions or information on the unknown dose–response function.

2.2 Spline-based point estimator of the target dose

Two components are needed to construct an estimator of the target dose. The first is an unbiased estimator of the response of the active control, which can be derived by using standard maximum likelihood (ML) theory (Pawitan 22) and can be written as the mean of the vector , that is, The more complex second part is to retrieve useful information regarding the unknown dose–response function between the investigated dose levels. Therefore, a smooth and continuous interpolation function without additional model assumptions is required. To construct a function with these attributes, a natural cubic spline interpolation based on the dose levels d1

The natural cubic spline function can be defined as with iteratively derived base polynomials and To derive the k + 2 base polynomials B, for q = 3, a knot vector of k + 6 data points is needed and can be defined by using multiples of the k dose levels or by using the mean difference between all dose levels . For the second approach, the knot vector can be defined as with d=d1+(i − 1)h, i =− q − 1,…,0 as well as d=d+jh, j = 1,…,q.

The natural cubic spline is generated by a sum of locally defined base polynomials (the B-spline basis) (de Boor 16) and connects each adjacent pair of dose levels d,d,i = 1,…,k − 1 with a cubic polynomial and has continuous first and second derivatives that guaranties a smooth curve shape. Further are the k conditions as well as the two natural spline conditions s′′(d1) = s′′(d) = 0 that are used to solve the unique linear system of k + 2 equations to estimate the unknown parameters c. This kind of spline interpolation is implemented in standard statistic software but can also be computed directly. When the parameters are derived, the point estimator of the target dose is derived numerically by solving the equation This is a point estimator of the target dose, which can be derived without any assumption on the shape of the dose–response function (except the assumptions made in section 1 which justify the use of the cubic spline).

2.3 Spline-based bootstrapped confidence interval for the target dose

In the previous section, a spline-based estimator of the target dose was presented. In the following section, the focus will lie on the construction of a CI for the target dose estimator. As described in section 1, no assumptions were made regarding the shape of the dose–response curve, and therefore, the distribution of the target dose estimator is unknown. To solve this problem, a bootstrap approach as well as the assumptions on the error terms in section 1 will be used. Even if the dose–response function is unknown, the mean values are unbiased estimators of f(d)i = 1,…,k, and the distributions of the mean values are known to be and with the active control (ac) labeled as the k + 1-th dose included in the study with n=n and σ=σ. The variances of the error terms can be estimated by Under these conditions, a bootstrap can be used to construct a CI. The residuals of the dose means will be bootstrapped using the variance estimator , and a new set of dose means , i = 1,…,k + 1 can be generated by (see, e.g., 17) with bootstrapped residuals and residual variance . As a second step, the interpolation method presented in section 2 will be used to compute target dose estimator , b = 1,…,n based on the bootstrapped mean values. All estimators of the target dose are saved in the vector , whose elements can then be ordered by size. In the last step, the α/2 and (1 − α/2) quantiles of are determined in order to construct the bootstrap CI. Then the CI based on the bootstrapped mean values and the cubic spline interpolation can be written as

3 Simulation studies

In this section, the properties of the presented methods from section 2 will be compared with standard parametric regression as well as nonparametric methods and be evaluated in extensive simulation studies motivated by the study presented in section 1. Therefore, we consider a five-arm dose-finding study (k = 4) with placebo (d1=0) as minimum dose and d4=1.8 as maximum dose ( = (d1,…,d4)′=(0,0.6,1.2,1.8)′). Further, different expected values of the active control are included (μ = 0.8,1.3) to enable the investigation of various target doses d*. The samples sizes per dose will be set to n=n, i = 1,…,k and the sample size of the active control to n=r·n with (r = 1). For all simulation scenarios, the number of replications is n=10,000, and the number of bootstrap samples per simulation run is n=5000, if not stated otherwise. Three different classes of dose–response functions are used in the simulation studies, the linear model f(d) = θ0+θ1d, the three parameter Emax model, and the four parameter logistic model (see Table I), which is sometimes also referred to as sigmoid Emax model 23. In the following, the term logistic model will be used to describe this response function.

Table I

Simulation scenarios for linear, Emax, and logistic response functions with θ0=−0.4, σ = 1.8, and sample sizes n = 20,25,…,40,50,…,100 per group.

Dose–response function	μ	θ1	θ2	θ3	d*
Linear	0.8	1.25	−	−	0.96
	1.3	1.25	−	−	1.36
Emax	0.8	2.675	0.4523	−	0.368
	1.3	2.675	0.4523	−	0.789
Logistic	0.8	2.675	0.9	3	0.840
	1.3	2.675	0.9	3	1.083

The Emax dose–response function can be described as with intercept θ0, slope θ1, and the dose, which leads to an effect of θ1/2, the ED50θ2. The logistic response function can be written as with intercept θ0, slope θ1, the ED50θ2, and the sigmoid shape parameter θ3. For d = 0, the response function is defined as f(0) = θ0, and the Emax model is a special case with θ3=1.

Two parametric nonlinear regressions (A, B) and a linear spline approach (C) presented by Dilleen et al. 15 are studied to compare the results of the cubic splines approach in the different dose–response scenarios. For both parametric regressions, the target dose estimator is derived by solving the equation , and the normal approximation by the Δ-method (Cramér's theorem, Ferguson 18) is applied to derive a CI for the target dose. (A) The first regression fits an Emax function and uses ML theory for estimating the response function and the target dose. (B) The second regression fits a logistic function and computes the target dose estimator as described for the first regression. To avoid fitting problems that can occur if a more flexible logistic regression model is fitted onto an Emax dose–response function, a three-step fitting procedure is used as described in Kirby et al. 24 and Jones et al. 25. If the logistic model cannot be fitted, an Emax model is fitted, and if even this was not possible, a linear model is used. (C) Furthermore, a flexible method based on linear spline interpolation presented by Dilleen et al. 15 was studied, which uses linear spline interpolation between the mean values of the dose levels for estimating the target dose. Dilleen et al. 15 assumed monotone dose–response functions and used isotonic regression to guarantee the monotonicity before nonparametric bootstrap methods were used to generate the CIs. This assumption of monotonicity is not made here, because the target dose was defined as the smallest dose that leads to the same efficacy as the active control and therefore is unique if it exists. To facilitate comparison of this approach versus the cubic spline approach of section 2, the bootstrap method described in section 3 is used instead of a nonparametric bootstrap approach, which was proposed by Dilleen et al. 15, to construct a CI for the linear spline approach. The three methods (A, B, and C) as well as the cubic splines of section 2 will be used in every simulation scenario, and the results will be presented in the next section.

3.1 Bias of the target dose estimator

In this section, the bias of the target dose estimators of the different methods will be investigated. The bias is defined as where is an estimator of d*. The ML estimators for the parametric regressions cannot be expressed in closed form and have to be calculated numerically. This problem is inherited by the bias in these cases, and therefore, the bias can be approximated using second-order Taylor approximations (see, for example, 26,27). Further, the bias of the ML-estimators in the case of normal error terms for N→∞ goes to zero. For the linear and the cubic spline, the bias is independent of the sample size, especially the bias does not decline to zero for N→∞. Therefore, we will focus on the investigation of the bias of the spline approaches. For the linear and the cubic spline interpolation, the bias cannot be given analytically but can be calculated numerically. The interpolations use the mean values of the different dose levels and the active control with to derive the target dose estimator. The mean values are unbiased estimators of f(d),i = 1,…,k, and therefore, the expected value of the spline interpolation can be derived by using the values of the expected response function f(d), i = 1…,k in the interpolation Hence, the bias is defined as The bias for the interpolations can be derived for every expected value μ of the active control to investigate the behavior of the bias over the dose range. The results of the numerically derived bias of the linear and the cubic spline as well as the corresponding spline interpolations for the Emax model are shown in Figure 2. It can be seen that the bias for the used active control responses (μ = 0.8,1.3) is 0.066 and −0.056 for the cubic spline and 0.104 and 0.063 for the linear spline, respectively. The maximum absolute bias over all reasonable μ∈[f(0),f(1.8)] is 0.087 for the cubic and 0.125 for the linear spline.

Figure 2

Linear and cubic spline interpolation of the Emax model fitted to the means reported in 21 and the bias of for the linear and cubic splines.

3.2 Coverage probability

The coverage probability of the 95% CI of d* using the different methods will be investigated through simulations in this section. The simulation specifications presented earlier as well as the simulation scenarios from Table I will be used.

To determine the performance of the different methods for the nominal coverage probability of 95%, the dotted lines in Figures 3 and 5 show the 99% interval of the simulation error and were derived as . This is roughly ±0.005, which has been considered as a practical irrelevant deviation from the nominal level 28. Simulated coverage probabilities above this interval are considered as conservative and values below as liberal.

Figure 3

Coverage probability of the cubic spline, the linear spline, the Emax, and logistic regression for different μ and various sample sizes n per group for a linear, an Emax and a logistic response function. The dotted lines indicate the simulation error with 99% probability. For all simulation scenarios, the number of replications is n=10,000 and the number of bootstrap simulations per simulation run is n=5000.

Figure 5

Cubic spline optimal: coverage probability of the cubic spline, the linear spline, the Emax, and logistic Δ-method for different μ and various sample sizes n per group for the Emax response function. The dotted lines indicate the simulation error with 99% probability. For all simulation scenarios, the number of replications is n=10,000 and the number of bootstrap simulations per simulation run is n=5000.

From Figure 3, it can be seen that for μ = 1.3, the cubic spline approach as well as the linear spline perform well for all sample sizes considered in the linear, the Emax, and the logistic response scenarios. Furthermore, it can be seen that the parametric Emax regression is working quite well for the linear and the Emax response function over all considered sample sizes, but in the cases of a logistic response function, it is not able to fit the underlying function and gets extremely liberal down to 88% for large sample sizes. The logistic regression is slightly liberal for all simulation scenarios with a linear or an Emax response function and performs well if the underlying response function is logistic. For an active control of μ = 0.8, all four methods are working acceptably for the linear and the logistic response function, but if an Emax response function was simulated, no method is able to hold the 95% coverage level. The Emax regression is liberal for small sample sizes and is getting better with increasing n, the linear spline and the logistic regression are getting extremely liberal and the effect is getting worse for large sample sizes. The same trend can be seen for the cubic spline but with a much smaller loss in coverage probability. For the spline approaches, this can be explained by the substantial bias in this scenario (see Figure 2) and the low number of probably unfavorable chosen dose levels. As described in Section 3, a three-step fitting procedure was used if the logistic model could not be fitted. In the investigated Emax scenarios, 30–38% of the logistic models could not be fitted over all sample sizes. In the Supporting Information, dose-finding studies with k = 5,6, and 8 dose levels have been investigated for the simulation scenarios described in Table I. Overall, the performance of the spline approaches improved with increasing number dose levels, whereas the parametric regression approaches do not improve with increasing number of dose levels and lead to liberal coverage probabilities, especially for small sample sizes.

In addition to the coverage probability, the median CI length has been investigated in the simulations studies. Under the considered scenarios, half open CIs can occur for the linear and the cubic spline interpolation because the spline approaches presented here are not designed to extrapolate beyond the dose range. These situations appear if, for example, the target dose estimator is close to the limits of the dose range and the sample sizes are small. The parametric regressions did not have this disadvantage, because they are able to extrapolate beyond the dose range. In the scenarios without half open CIs, all methods perform very similar for all considered response functions (results not shown).

3.3 The cubic spline with penalty term and cross-validation

The presented spline-based estimation of the target dose and the combination with the bootstrapped CI can be used with every kind of fitting approach as long as no other assumptions are needed. Therefore, the potential benefit of a more complex smoothing cubic penalty spline as proposed for dose-finding studies in 14 is investigated in the following. Instead of interpolating the mean responses at the dose levels, a cubic spline in combination with a penalty term is used and the following equation has to be optimized where α is the so-called smoothing parameter and the second derivative of the cubic spline function. The α has to be prespecified or can, for example, be derived via cross-validation, which was used in 14. To compare this more complex spline approach with the presented cubic spline interpolation, a small simulation study is conducted, which uses the presented simulation scenarios of the Emax response function in Table I for k = 4,8 dose levels and two different sample sizes n = 25,50 per group with n=2000 and n=1000. The designs with k = 8 dose levels were included in the simulation study to address the potential advantage of the penalty spline, which is the smoothing of noisy data. Additionally, to the cubic spline interpolation and the selection of the penalty term by cross-validation, two prespecified penalty terms for medium and high penalty (in SAS 9.3 log10(k·α) =− 4,−2) are included in the simulations. To derive the cubic penalty spline, the standard procedure TPSPLINE in SAS 9.3 was applied, which is able to use predefined α as well as to derive α by cross-validation. The results are shown in Table II. It can be seen that there is no substantial benefit of the penalty spline over the cubic interpolation in the investigated scenarios. Further, it can be seen that choosing the penalty term α by cross-validation leads to a loss in coverage probability in the case of k = 4 dose levels, which can be explained by the low number of dose levels. In these cases, it is more useful to predefine the penalty term or use the presented spline interpolation. For designs with k = 8 dose levels, the spline interpolation as well as penalty splines with selected penalty terms and with cross-validation perform very well. In terms of target dose estimation and coverage probability, no real benefit of the penalty spline over the spline interpolation was found. Additionally, in designs with a low number of dose levels, the performance of the penalty spline is fairly poor in comparison with the spline interpolation approach. This was one of the reasons to use the simpler and more stable cubic spline interpolation. If the interest of the study also lies on the identification of the dose–response shape, the penalty spline can compute much smoother curves especially for higher numbers of dose levels.

Table II

The mean estimated target dose, coverage probability, and median confidence interval length for different numbers of dose levels (k), different expected responses of the active control (μ), and different sample sizes (n) per group are reported for the cubic spline interpolation (penalty = no) and the smoothing spline with various penalty terms (penalty = cross-validation (CV), High, Medium) in the Emax response scenario (Table I) with n=2000 and n=1000.

k	μ	d*	n	Penalty	Mean	Coverage probability	Median interval length
				CV	0.5814	0.8875	1.2090
			25	High	0.5301	0.9285	1.1870
				Medium	0.5089	0.9305	1.1775
	0.8	0.3680		No	0.4931	0.9350	1.1680
				CV	0.5267	0.8750	0.8510
			50	High	0.5023	0.8825	0.7360
				Medium	0.4651	0.9175	0.6975
				No	0.4578	0.9260	0.7000
4				CV	0.9722	0.9490	1.3590
			25	High	0.9302	0.9425	1.4040
				Medium	0.8814	0.9465	1.4365
	1.3	0.7886		No	0.8753	0.9410	1.4350
				CV	0.9197	0.9505	1.2100
			50	High	0.9032	0.9570	1.2525
				Medium	0.8530	0.9465	1.3220
				No	0.8577	0.9435	1.3150
				CV	0.4668	0.9560	1.1405
			25	High	0.4723	0.9295	1.0960
				Medium	0.4205	0.9485	0.9170
	0.8	0.3680		No	0.4278	0.9555	0.9040
				CV	0.4388	0.9520	0.7900
			50	High	0.4546	0.9335	0.7030
				Medium	0.3962	0.9560	0.7070
				No	0.4030	0.9455	0.6960
8				CV	0.9134	0.9525	1.4495
			25	High	0.9074	0.9490	1.4090
				Medium	0.7960	0.9385	1.4010
	1.3	0.7886		No	0.7734	0.9420	1.3860
				CV	0.8503	0.9705	1.2790
			50	High	0.8751	0.9435	1.2695
				Medium	0.8069	0.9545	1.2325
				No	0.7976	0.9455	1.2190

4 Designs for the spline-based approach: minimizing the bias of the target dose estimator

Now, we investigate the design of a k dose level dose-finding study with active control. The minimum and the maximum dose levels (d1,d) are derived from the information by tolerability phase I studies or determined in some other way by the clinical trial team. Also, the ratio r between dose sample sizes n=n, for i = 1,…,k per group and active control sample size n=r·n is prespecified. Even if the spline approaches are very flexible, some information or assumptions about the underlying dose–response function are needed in the planning phase. So the optimal designs are depending on the class of the dose–response function (e.g., linear, Emax or logistic model) and will be derived in three different ways as presented in the succeeding text. The aim of all procedures is to find a design, which is optimal for an unknown response of the active control so that the absolute bias of the spline procedures is small overall.

4.1 Minimum maximum global bias

Using the Taylor series approximation it was shown in Hall 29, Hall and Meyer 30 and de Boor 16 that the spline interpolations have the following global upper error bounds. For the linear spline s1, the upper error bound is given by and for the cubic spline s3 the upper error bound is with constants , , and as well as . If we want to minimize these upper bounds, h must be minimized and this can be achieved by choosing the dose levels equidistant between d1 and d so that h=d−d for i = 1,…,k − 1. With this approach, we minimize the general upper bounds for the spline. However, this approach is mainly used to derive asymptotic results regarding the number of dose levels going to infinity (k→∞) (see de Boor 16). This is not the goal of the designs we have in mind because we consider asymptotic properties with n→∞ or N→∞ and a fixed number of dose levels. Because this optimization does not depend on the underlying dose–response function, it can be done with no extra information.

4.2 Equally spaced responses

In this approach, the dose levels will be selected in a way that the dose–response difference between two adjacent dose levels is the same for all dose levels With this condition, more information is available on the parts of the response function with larger changes in the response. For the linear, Emax, and logistic model, it can be shown that this approach is independent of the intercept θ0 and the slope θ1 (see Appendix A). This means that for the linear model, an equidistant allocation of the dose levels on the dose range would be optimal and that no assumptions have to be made to generate the optimal design. This result leads also to the same optimal design for the linear model as the global bias approach in the previous section. Furthermore, for the Emax model, only the ED50(θ2) must be specified to derive the optimal design for the linear as well as the cubic spline. In the logistic model, additionally to the ED50, the degree of sigmoidal shape (θ3) must be quantified.

4.3 Equally bound bias

This approach uses the global upper error bound of the bias in Section 4.1 on every interval between two adjacent dose levels simultaneously for all dose levels. Therefore, we define the optimization criteria for the linear spline s1 as and for the cubic spline s3 as Using this approach, the bias is equally bound on every interval between two adjacent dose levels ensuring a good interpolation over the investigated dose range. As in the approaches presented earlier, it can be shown that this optimization procedure does not depend on the intercept θ0 and the slope θ1 for linear, Emax, and logistic models (see Appendix A). The construction of the optimal designs can be done in the same way as described in Section 4.2. It is important to note that this method as well as the approach described in the previous section need numerical optimization to find the corresponding optimal design. This can be done in any standard statistical software (e.g., SAS, R) using optimization routines for solving root-finding problems. To use this kind of standard routines for the presented methods, all k − 1 conditions can be arranged in a vector and must be multiplied by the projection matrix . Then the optimization is a root-finding problem with the same result, and the standard routines can be applied.

4.3.1 The minimum–maximum global absolute bias

To determine for the desired spline interpolation, which of the methods (equidistant on the dose scale (m = 1), equidistant on the response scale (m = 2), equally bound bias (m = 3)) is superior to the others, it is necessary to define a criterion that can assess the differences between the methods. Such a criterion is the minimum–maximum global absolute bias of the spline interpolation, which is given by with Bias the bias of method m.

Through this criterion, the method m with the smallest maximum absolute bias over all reasonable expected values of the active control μ(c,c∈[0,(f(d) − f(d1))/2]) is chosen as the optimal design. Following this line, a design can be found that is optimal for all reasonable μ and has the smallest bias of the presented methods. If the study is placebo controlled (d1=0), then the already established active comparator should have a higher response than placebo, and c>0 would be a reasonable choice. In our example, c and c are set to 0.

In Figure 4, the bias-based optimal design for the cubic spline is presented for the Emax model used in the simulation studies, with θ2=0.4523 to enable comparison of the bias in the original design (see Figure 2). It can be seen that in the case of four dose levels, the optimal design for the cubic spline leads to a larger bias for certain expected values in the linear spline interpolation than in the equidistant case. In the simulation scenarios μ = 0.8,1.3, the bias in the optimal design is −0.017 and 0.020 for the cubic spline and 0.037 and 0.223 for the linear spline, respectively. The maximal absolute bias is then 0.023 for the cubic spline and 0.267 for the linear spline. In such situations, we would not recommend the use of linear splines in optimal designs which were planned for the use of cubic splines. If more than four dose levels are investigated, the bias decreases further for both methods, and the optimal designs for the linear and the cubic spline become more similar.

Figure 4

Planning phase: linear and cubic spline interpolation of the Emax model fitted to the means reported in 21 and the bias of for the linear and the cubic spline in the bias optimal design of the cubic spline.

4.3.2 Bias minimal design: simulation study

To investigate the influence of the reduced bias approaches on the coverage probability, new simulations were performed in the bias minimal Emax designs for the simulation scenarios of the Emax response function with θ2=0.4523 (Table I). The Emax model with the new design =(d1,d2,…,d4)′=(0,0.1646,0.516,1.8)′ was simulated with n=10,000 and n=5000 replications for the sample sizes n = 20,25,30,40,…,100 per group. The results are presented in Figure 5.

It can be seen that for μ = 0.8, the coverage probability of the cubic spline approach stays in comparison with Figure 3 within the 99% interval of the simulation error for all sample sizes. For μ = 1.3, the coverage probability becomes a little bit conservative. The bias minimizing design leads not only to a small bias but also to a much better coverage probability in the investigated scenarios. However, this is only the case for the cubic spline. The performance of the linear spline does not improve because the bias of the linear spline is still very large in the optimal design for a cubic spline. To be fair, the optimal design for the linear spline is calculated as well and leads in this example to the equidistant allocation, which was already used in the example and the simulation study. Therefore, no extra simulation study is executed and for results, it is referred to Section 3. It is also interesting to see that both parametric models using the Δ-method are getting liberal in the simulation scenarios and are therefore also sensitive to dose selection.

5 Conclusions

Extending the previous work by Dilleen et al. 15, cubic spline interpolation was proposed for the purpose of estimating the target dose and constructing bootstrapped CIs. To investigate the performance of this approach, simulation studies were performed. We could show that the cubic spline approach has a reasonable coverage probability for the linear, the Emax, and the logistic response function in most of the scenarios considered. For the more biased case of an Emax response with expectation of the active control μ = 0.8, we saw that the cubic spline is slightly liberal and this tendency increases with larger sample sizes n per group. The performance of the linear spline approach is quite similar to the cubic spline interpolation, but in the described Emax scenarios with μ = 0.8, the coverage probability becomes extremely liberal.

A key factor in explaining these problems is the large bias of the linear and the cubic spline, which becomes increasingly dominant with decreasing standard errors for larger sample sizes n per group. Additionally, the resulting CI length for μ = 0.8 is much shorter than for the case of μ = 1.3. Even if the absolute bias is comparable in the two scenarios, the coverage probability is more sensitive to the underlying bias if the interval length is shorter. The difference of the interval length can be seen in Table II. The cubic spline interpolation is represented by (penalty = no), for n = 50 and k = 4, the median interval length is 0.7 for μ = 0.8 and 1.31 for μ = 0.8, which is nearly the entire dose range. The Emax regression works fine if a linear or an Emax response function was simulated and becomes very liberal if the underlying function for the simulation was a logistic function. The coverage probability of the logistic regression is liberal for nearly all simulation scenarios, except if the simulated response function was a logistic function. Furthermore, in 30–38% of the simulation runs, it was not possible to fit a logistic regression when an Emax response function was simulated. Instead of using the three-step fitting procedure, parameter bounds could be used in the regression procedures to avoid the numerical problems (see, for example, 31). Additionally to the designs with k = 4 dose levels, designs with k = 5,6, and 8 dose levels were simulated, which are reported in the Supporting Information. Because of the construction of the spline approaches, an extrapolation beyond the dose range was not possible and one-sided intervals occurred. In these cases, the median length was infinity, and therefore, a comparison was not sensible. On the other hand, if the median length was finite for all approaches, no large differences could be seen between the four approaches. To justify the use of the cubic spline interpolation instead of a cubic smoothing spline with cross-validation, a small simulation study was conducted. It could be shown that the coverage probability of the smoothing spline is quite liberal if only four dose levels were investigated. Overall, the cubic spline interpolation performed much better in the simulated scenarios. Therefore, in the cases of small numbers of dose levels, the cubic spline interpolation was preferred.

Different bias reducing procedures were described and combined to find an optimal design by the criteria of the minimal maximum absolute bias on the desired active control (or response) range. The planning of such studies is relatively simple because only very few parameters must be prespecified to construct the optimal design. For instance, for an Emax response function, only the ED50=θ2 and for a logistic function, only the ED50=θ2 and the sigmoidal shape parameter θ3 must be prespecified. In the simulation study for the optimal design (Section 4.3.2), the true θ2 value of the presented Emax response function was used to derive the new dose allocation. We could show that the bias reducing design leads to a considerably better control of the coverage probability. On the other hand, the coverage probability of parametric approaches is much too low especially for smaller and moderate sample sizes n = 20,…,70 per group. The performances of the parametric methods depend heavily on the investigated dose levels and could lead to extreme liberal results in the investigated scenarios. With an increased number of dose levels, the splines become more stable and less dependent on the assumed parameters. We recommend the use of the spline-based methods only in studies with at least four dose levels in situations where some knowledge of the true parameters exists or five dose levels to achieve a reasonable performance. For the calculation of the optimal spline design, some of the model parameters of the response function are needed and usually unknown during the planning phase. Therefore, the extension of the optimal spline design under model uncertainty would be very interesting. A Min–Max approach could be used to derive the design not only for one fix dose response function but also for multiple parameter constellations. The resulting dose allocation would minimize the maximum bias over all considered model variations. In addition, the use of optimal designs for cubic splines for the planning of studies that should be evaluated using parametric approaches such as the Δ-method would be of further interest, especially in the context of comparing the spline-based optimal designs with the optimal designs used in parametric approaches, which are based on other local optimality criteria (see, for example, Dette et al. 13,32).