Optimal designs for model averaging in non-nested models

TL;DR

This paper develops optimal experimental designs for model averaging in non-nested, misspecified models, showing that Bayesian optimal designs significantly enhance estimation accuracy.

Contribution

It introduces Bayesian optimal designs for model averaging in non-nested models with fixed weights, improving estimation accuracy over traditional methods.

Findings

Bayesian optimal designs reduce mean squared error of estimates.

Designs improve accuracy for both model averaging and model selection.

Optimal designs outperform non-optimized designs in simulations.

Abstract

In this paper we construct optimal designs for frequentist model averaging estimation. We derive the asymptotic distribution of the model averaging estimate with fixed weights in the case where the competing models are non-nested and none of these models is correctly specified. A Bayesian optimal design minimizes an expectation of the asymptotic mean squared error of the model averaging estimate calculated with respect to a suitable prior distribution. We demonstrate that Bayesian optimal designs can improve the accuracy of model averaging substantially. Moreover, the derived designs also improve the accuracy of estimation in a model selected by model selection and model averaging estimates with random weights.

Figures

Tables

Table 1. Table 1: Models and parameters used for the simulation study.

Model	Mean function ${\eta }_s$	Parameter specifications
Log-Linear ($f_1$)	${\eta }_1(x_i,{\vartheta }_1)={\vartheta }_{11}+{\vartheta }_{12}\log \left(x_i+{\vartheta }_{13}\right)$	${\vartheta }_1={(0,0.0797,1)}^{\top }$
Emax ($f_2$)	${\eta }_2(x_i,{\vartheta }_2)={\vartheta }_{21}+\frac{{\vartheta }_{22}x}{{\vartheta }_{23}+x}$	${\vartheta }_2={(0,0.467,25)}^{\top }$
Exponential ($f_3$)	${\eta }_3(x_i,{\vartheta }_3)={\vartheta }_{31}+{\vartheta }_{32}\exp \left(x_i/{\vartheta }_{33}\right)$	${\vartheta }_3={(-0.08265,0.08265,85)}^{\top }$
Quadratic ($f_4$)	${\eta }_4(x_i,{\vartheta }_4)={\vartheta }_{41}+{\vartheta }_{42}x+{\vartheta }_{43}{x}^2$	${\vartheta }_4={(0,0.00533,-0.00002)}^{\top }$

Table 2. Table 2: Simulated mean squared error of different estimates of the ED 0.4 subscript ED 0.4 \text{ED}_{0.4} . The set of candidate models is 𝒮 1 = { f 1 , f 2 , f 4 } subscript 𝒮 1 subscript 𝑓 1 subscript 𝑓 2 subscript 𝑓 4 {\cal S}_{1}=\{f_{1},f_{2},f_{4}\} . Left column: model averaging with uniform weights; middle column: model averaging with smooth AIC-weights; right column: estimation after model selection.

model	sample size	uniform weights	smooth AIC-weights	model selection
	$n=50$	437.045	498.323	758.978
$f_1$	$n=100$	223.291	218.99	285.062
	$n=250$	111.973	82.713	78.371
	$n=50$	286.638	329.904	515.32
$f_2$	$n=100$	189.785	203.796	251.836
	$n=250$	62.792	64.854	66.54
	$n=50$	276.037	361.101	669.873
$f_4$	$n=100$	190.662	244.558	391.443
	$n=250$	92.653	109.852	139.859
	$n=50$	1503.903	1372.31	1381.033
$f_3$	$n=100$	1109.622	856.484	729.912
	$n=250$	864.163	398.144	255.604

Table 3. Table 3: Simulated mean squared error of different estimates of the ED 0.4 subscript ED 0.4 \text{ED}_{0.4} . The set of candidate models is 𝒮 2 = { f 1 , f 2 , f 3 } subscript 𝒮 2 subscript 𝑓 1 subscript 𝑓 2 subscript 𝑓 3 {\cal S}_{2}=\{f_{1},f_{2},f_{3}\} . Left column: model averaging with uniform weights; middle column: model averaging with smooth AIC-weights; right column: estimation after model selection.

		estimation method
model	sample size	uniform weights	smooth AIC-weights	model selection
	$n=50$	834.295	553.427	776.311
$f_1$	$n=100$	712.404	340.254	353.707
	$n=250$	524.518	48.587	38.591
	$n=50$	640.706	505.054	669.285
$f_2$	$n=100$	517.963	267.967	286.272
	$n=250$	394.536	65.805	53.424
	$n=50$	1076.154	1141.476	1427.441
$f_3$	$n=100$	871.362	766.140	802.763
	$n=250$	802.196	480.641	399.839
	$n=50$	288.091	486.501	852.377
$f_4$	$n=100$	208.628	298.315	419.651
	$n=250$	162.689	138.331	142.673

Table 4. Table 4: Simulated mean squared errors of different estimates of the ED 0.4 subscript ED 0.4 \text{ED}_{0.4} for different experimental designs. The set of candidate models is 𝒮 1 = { f 1 , f 2 , f 4 } subscript 𝒮 1 subscript 𝑓 1 subscript 𝑓 2 subscript 𝑓 4 \mathcal{S}_{1}=\{f_{1},f_{2},f_{4}\} . Left column: model averaging estimate with uniform weights; middle column: model averaging estimate with smooth AIC-weights; right column: estimate after model selection.

model	design	uniform weights	smooth AIC-weights	model selection
	(4.2)	177.472	165.981	173.548
$f_1$	(4.3)	223.291	218.990	285.062
	(3.25)	185.251	184.77	340.698
	(4.2)	142.085	153.745	170.059
$f_2$	(4.3)	189.785	203.796	251.836
	(3.25)	501.814	501.394	1162.654
	(4.2)	160.039	195.116	299.365
$f_4$	(4.3)	190.662	244.558	391.443
	(3.25)	404.716	427.548	1396.051
	(4.2)	1058.655	766.408	606.752
$f_3$	(4.3)	1109.622	856.484	729.912
	(3.25)	3184.11	3413.566	4102.964

Table 5. Table 5: Simulated mean squared errors of different estimates of the ED 0.4 subscript ED 0.4 \text{ED}_{0.4} for different experimental designs. The set of candidate models is 𝒮 2 = { f 1 , f 2 , f 3 } subscript 𝒮 2 subscript 𝑓 1 subscript 𝑓 2 subscript 𝑓 3 \mathcal{S}_{2}=\{f_{1},f_{2},f_{3}\} . Left column: model averaging estimate with uniform weights; middle column: model averaging estimate with smooth AIC-weights; right column: estimate after model selection.

		estimation method
model	design	uniform weights	smooth AIC-weights	model selection
	(4.4)	654.914	279.257	274.016
$f_1$	(4.3)	712.404	340.254	353.707
	(3.25)	770.705	410.715	413.676
	(4.4)	546.105	283.757	250.719
$f_2$	(4.3)	517.963	267.967	286.272
	(3.25)	1098.323	962.257	1701.569
	(4.4)	910.372	742.507	699.612
$f_3$	(4.3)	871.362	766.140	802.763
	(3.25)	1505.693	1774.895	2592.261
	(4.4)	159.899	278.409	347.187
$f_4$	(4.3)	208.628	298.315	419.651
	(3.25)	522.652	610.198	1907.066

Equations

$$Y_{ij}={\eta }_s(x_i,{\vartheta }_s)+{\epsilon }_{ij},i=1,\dots ,k,j=1,\dots ,n_i,s=1,\dots ,r,$$

$${\hat{\theta }}_{n,s}=\text{arg}\underset{{\theta }_s\in {\Theta }_s}{\max }{\ell }_{n,s}({\theta }_s\mid y_n)$$

$${\ell }_{n,s}({\theta }_s\mid y_n)=\frac{1}{n}\sum _{i=1}^k\sum _{j=1}^{n_i}\log f_s(y_{ij}\mid x_i,{\theta }_s)$$

$$\mu_{s}(\theta_{s})=\inf\left\{x\in{\cal X}~{}\Big{|}~{}\tfrac{\eta_{s}(x,\vartheta_{s})-\eta_{s}(a,\vartheta_{s})}{\eta_{s}(b,\vartheta_{s})-\eta_{s}(a,\vartheta_{s})}\geq\alpha\right\}~{}.$$

$$\text{AIC}(f_s\mid y_n)=2{\ell }_{n,s}({\hat{\theta }}_{n,s}\mid y_n)-2p_s,$$

$${\hat{\mu }}_{\text{mav}}=\sum _{s=1}^r{w}_s{\mu }_s({\hat{\theta }}_{n,s}).$$

$$w_s^{\text{smAIC}}(y_n)=\frac{\exp \{\frac{1}{2}\text{AIC}(f_s\mid y_n)\}}{{\sum }_{s=1}^r\exp \{\frac{1}{2}\text{AIC}(f_s\mid y_n)\}}.$$

$${\mathcal{S}}_1=\{f_1,f_2,f_4\}$$

$${\mathcal{S}}_2=\{f_1,f_2,f_3\}$$

$$\frac{1}{6}\sum _{i=1}^6\int f_4(y\mid x_i,{\theta }_4)\log \left(\frac{f_4(y\mid x_i,{\theta }_4)}{f_s(y\mid x_i,{\theta }_s)}\right)dy,s=1,2.$$

$$y_{ij}^{(l)}={\eta }_s(x_i,{\vartheta }_s)+{\epsilon }_{ij}^{(l)},i=1,\dots ,k,j=1,\dots ,n_i,$$

$$f_1(y\mid x,{\theta }_1),\dots ,f_r(y\mid x,{\theta }_r)$$

$$\xi =\left\{x_1,\dots ,x_k;{\xi }_1,\dots ,{\xi }_k\right\},$$

$$\text{KL}(g:f_s\mid {\theta }_s,\xi )=\int g(y\mid x)\log \left(\frac{g(y\mid x)}{f_s(y\mid x,{\theta }_s)}\right)dyd\xi (x),$$

$${\theta }_{s,g}^{\ast }(\xi )=\text{arg}\underset{{\theta }_s\in {\Theta }_s}{\min }\text{KL}(g:f_s\mid {\theta }_s,\xi )$$

$$A_s({\theta }_s,\xi )$$

$$B_{st}({\theta }_s,{\theta }_t,\xi )$$

$$\sqrt{n}({\hat{\theta }}_{n,s}-{\theta }_s^{\ast }(\xi ))\stackrel{\mathcal{D}}{⟶}\mathcal{N}\left(0,A_s^{-1}({\theta }_s^{\ast }(\xi ))B_{ss}({\theta }_s^{\ast }(\xi ),{\theta }_s^{\ast }(\xi ))A_s^{-1}({\theta }_s^{\ast }(\xi ))\right),$$

$$A_s({\theta }_s^{\ast }(\xi ))=A_s({\theta }_s^{\ast }(\xi ),\xi ),B_{st}({\theta }_s^{\ast }(\xi ),{\theta }_t^{\ast }(\xi ))=B_{st}({\theta }_s^{\ast }(\xi ),{\theta }_t^{\ast }(\xi ),\xi )$$

$$\displaystyle\sqrt{n}\Big{(}\hat{\mu}_{\text{mav}}-\sum_{s=1}^{r}w_{s}\mu_{s}(\theta_{s}^{*}(\xi))\Big{)}\overset{\mathcal{D}}{\longrightarrow}\mathcal{N}\left(0,\sigma_{w}^{2}(\theta^{*}(\xi))\right),$$

$$\displaystyle\sigma_{w}^{2}(\theta^{*}(\xi))=\sum_{s,t=1}^{r}w_{s}w_{t}\Big{(}\tfrac{\partial\mu_{s}(\theta_{s}^{*}(\xi))}{\partial\theta_{s}}\Big{)}^{\top}A_{s}^{-1}(\theta_{s}^{*}(\xi))B_{st}\left(\theta_{s}^{*}(\xi),\theta_{t}^{*}(\xi)\right)A_{t}^{-1}\left(\theta_{t}^{*}(\xi)\right)\tfrac{\partial\mu_{t}(\theta_{t}^{*}(\xi))}{\partial\theta_{t}}.$$

$$\displaystyle\Phi_{\text{mav}}(\xi,g)=\tfrac{1}{n}\sigma_{w}^{2}(\theta^{*}(\xi))+\Big{(}\sum_{s=1}^{r}w_{s}\mu_{s}(\theta_{s}^{*}(\xi))-\mu_{\text{true}}\Big{)}^{2}\approx\text{MSE}(\hat{\mu}_{\text{mav}}),$$

$${\Phi }_{\text{mav}}^{\pi }(\xi )={\int }_{\mathcal{G}}{\Phi }_{\text{mav}}(\xi ,g)d\pi (g).$$

$$d\pi (g,\theta )=d{\pi }_1(\theta \mid g)d{\pi }_2(g),$$

$${\Phi }_{\text{mav}}^{\pi }(\xi )={\int }_{\mathcal{S}}{\int }_{{\mathcal{F}}_g}{\Phi }_{\text{mav}}(\xi ,g)d{\pi }_1(\theta \mid g)d{\pi }_2(g),$$

$$d_{\pi }(x,{\xi }^{\ast })=$$

$${\theta }_{s,g}^{\prime }({\xi }^{\ast },x)$$

$$\int g(y\mid x)\frac{\partial }{\partial {\theta }_s}\log f_s(y\mid x,{\theta }_{s,g}^{\ast }({\xi }^{\ast }))dy$$

$${\sigma }_g^{\prime }({\xi }^{\ast },x)$$

$$\displaystyle-\Big{(}\tfrac{\partial\mu_{s}(\theta_{s,g}^{*}(\xi^{*}))}{\partial\theta_{s}}\Big{)}^{\top}\big{(}(A_{s}^{*})^{-1}h_{s,g}^{\prime}(\xi^{*},x)(A_{s}^{*})^{-1}\big{)}B_{st}^{*}(A_{t}^{*})^{-1}\tfrac{\partial\mu_{t}(\theta_{t,g}^{*}(\xi^{*}))}{\partial\theta_{t}}$$

$$\displaystyle+\Big{(}\tfrac{\partial\mu_{s}(\theta_{s,g}^{*}(\xi^{*}))}{\partial\theta_{s}}\Big{)}^{\top}(A_{s}^{*})^{-1}h_{st,g}^{\prime}(\xi^{*},x)(A_{t}^{*})^{-1}\tfrac{\partial\mu_{t}(\theta_{t,g}^{*}(\xi^{*}))}{\partial\theta_{t}}$$

$$\displaystyle-\Big{(}\tfrac{\partial\mu_{s}(\theta_{s,g}^{*}(\xi^{*}))}{\partial\theta_{s}}\Big{)}^{\top}(A_{s}^{*})^{-1}B_{st}^{*}\left((A_{t}^{*})^{-1}h_{t,g}^{\prime}(\xi^{*},x)(A_{t}^{*})^{-1}\right)\tfrac{\partial\mu_{t}(\theta_{t,g}^{*}(\xi^{*}))}{\partial\theta_{t}}$$

$$\displaystyle+\Big{(}\tfrac{\partial\mu_{s}(\theta_{s,g}^{*}(\xi^{*}))}{\partial\theta_{s}}\Big{)}^{\top}(A_{s}^{*})^{-1}B_{st}^{*}(A_{t}^{*})^{-1}\tfrac{\partial^{2}\mu_{t}(\theta_{t,g}^{*}(\xi^{*}))}{\partial\theta_{t}\partial\theta_{t}^{\top}}\theta_{t,g}^{\prime}(\xi^{*},x)\Big{]},$$

$$h_{st,g}^{\prime }({\xi }^{\ast },x)$$

Full text

Keywords: Model selection, model averaging, model uncertainty, optimal design, Bayesian optimal design

1 Introduction

There exists an enormous amount of literature on selecting an adequate model from a set of candidate models for statistical analysis. Numerous model selection criteria have been developed for this purpose. These procedures are widely used in practice and have the advantage of delivering a single model from a class of competing models, which makes them very attractive for practitioners. Exemplarily, we mention Akaike’s information criterion (AIC), the Bayesian information criterion (BIC) and its extensions, Mallow’s $C_{p}$, the generalized cross-validation and the minimum description length (see the monographs of Burnham and Anderson, (2002), Konishi and Kitagawa, (2008) and Claeskens and Hjort, (2008) for more details). Different criteria have different properties, such as consistency, efficiency and parsimony (used in the sense of Claeskens and Hjort, (2008, Chapter 4)). Overall there seems to be no universally optimal model selection criterion and different criteria might be preferable in different situations depending on the particular application.

On the other hand, there exists a well known post-selection problem in this approach because model selection introduces an additional variance that is often ignored in statistical inference after model selection (see Pötscher, (1991) for one of the first contributions discussing this issue). This post-selection problem is inter alia attributable to the fact, that estimates after model selection behave like mixtures of potential estimates. For example, ignoring the model selection step (and thus the additional variability) may lead to confidence intervals with coverage probability smaller than the nominal value, see for example Chapter 7 in Claeskens and Hjort, (2008) for a mathematical treatment of this phenomenon.

An alternative to model selection is model averaging, where estimates of a target parameter are smoothed across several models, rather than restricting inference on a single selected model. This approach has been widely discussed in the Bayesian literature, where it is known as “Bayesian model averaging” (see the tutorial of Hoeting et al., (1999) among many others). For Bayesian model averaging prior probabilities have to be specified. This might not always be possible and therefore Hjort and Claeskens, (2003) also proposed a “frequentist model averaging”, where smoothing across several models is commonly based on information criteria. Kapetanios et al., (2008) demonstrated that the frequentist approach is a worthwhile alternative to Bayesian model averaging. Stock and Watson, (2003) observed that averaging predictions usually performs better than forecasting in a single model. Hong and Preston, (2012) substantiate these observations with theoretical findings for Bayesian model averaging if the competing models are “sufficiently close”. Further results pointing in this direction can be found in Raftery and Zheng, (2003), Schorning et al., (2016) and Buatois et al., (2018).

Independently of this discussion there exists a large amount of research how to optimally design experiments under model uncertainty (see Box and Hill, (1967); Atkinson and Fedorov, (1975) for early contributions). This work is motivated by the fact that an optimal design can improve the efficiency of the statistical analysis substantially if the postulated model assumptions are correct, but may be inefficient if the model is misspecified. Many authors suggested to choose the design for model discrimination such that the power of a test between competing regression models is maximized (see Ucinski and Bogacka, (2005); López-Fidalgo et al., (2007); Tommasi and López-Fidalgo, (2010) or Dette et al., (2015) for some more recent references). Other authors proposed to minimize an average of optimality criteria from different models to obtain an efficient design for all models under consideration (see Dette, (1990), Zen and Tsai, (2002); Tommasi, (2009) among many others).

Although model selection or averaging are commonly used tools for statistical inference under model uncertainty most of the literature on designing experiments under model uncertainty does not address the specific aspects of these methods directly. Optimal designs are usually constructed to maximize the power of a test for discriminating between competing models or to minimize a functional of the asymptotic variance of estimates in the different models. To the best of our knowledge Alhorn et al., (2019) is the first contribution, which addresses the specific challenges of designing experiments for model selection or model averaging. These authors constructed optimal designs minimizing the asymptotic mean squared error of the model averaging estimate and showed that optimal designs can yield a reduction of the mean squared error up to $45\%$. Moreover, they also showed that these designs improve the performance of estimates in models chosen by model selection criteria. However, their theory relies heavily on the assumption of nested models embedded in a framework of local alternatives as developed by Hjort and Claeskens, (2003).

The goal of the present contribution is the construction of optimal designs for model averaging in cases where the competing models are not nested (note that in this case local alternatives cannot be formulated). Moreover, in contrast to most of the literature, we also consider the situation where all competing models misspecify the data underlying truth. In order to derive an optimality criterion, which can be used for the determination of optimal designs in this context, we further develop the approach of Hjort and Claeskens, (2003) and derive an asymptotic theory for model averaging estimates for classes of competing models which are non-nested. Optimal designs are then constructed minimizing the asymptotic mean squared error of the model averaging estimate and it is demonstrated that these designs yield substantially more precise model averaging estimates. Moreover, these designs also improve the performance of estimates after model selection. Our work also contributes to the discussion of the superiority of model averaging over model selection. Most of the results presented in literature indicate that model averaging has some advantages over model selection in general. We demonstrate that conclusions of this type depend sensitively on the class of models under consideration. In particular we observe some advantages of estimation after model selection if the competing models are of rather different shape. Nevertheless, the optimal designs developed in this paper improve both estimation methods, where the improvement can be substantial in many cases.

The remaining part of this paper is organized as follows. The pros and cons of model averaging and model selection are briefly discussed in Section 2 where we introduce the basic methodology and investigate the impact of similarity of the candidate models on the performance of the different estimates. In Section 3 we develop asymptotic theory for model averaging estimation in the case where the models are non-nested and all competing models might misspecify the underlying truth. Based on these results we derive a criterion for the determination of optimal designs and establish a necessary condition for optimality. In Section 4 we study the performance of these designs by means of a simulation study. Finally, technical assumptions and proofs are given Section 6.

2 Model averaging versus model selection

In this section we introduce the basic terminology and also illustrate in a regression framework that the superiority of model averaging about estimation in a model chosen by model selection depends sensitively on the class of competing models.

2.1 Basic terminology

We consider data obtained at $k$ different experimental conditions, say $x_{1},\ldots,x_{k}$ chosen in a design space $\mathcal{X}$. At each experimental condition $x_{i}$ one observes $n_{i}$ responses, say $y_{i1},\ldots,y_{in_{i}}$ $(i=1,\ldots,k)$, and the total sample size is $n=\sum_{i=1}^{k}n_{i}$. We also assume that the responses $y_{i1},\ldots,y_{in_{i}}$ are realizations of random variables of the form

[TABLE]

where the regression function $\eta_{s}$ is a differentiable function with respect to the parameter $\vartheta_{s}$ and the random errors $\varepsilon_{ij}$ are independent normally distributed with mean 0 and common variance $\sigma^{2}$. Furthermore, the index $s$ in $\eta_{s}$ corresponds to different models (with parameters $\vartheta_{s}$) and we assume that there are $r$ competing regression functions $\eta_{1},\ldots,\eta_{r}$ under consideration.

Having $r$ different candidate models (differing by the regression functions $\eta_{s}$) a classical approach for estimating a parameter of interest, say $\mu$, is to calculate an information criterion for each model under consideration and estimate this parameter in the model optimizing this criterion. For this purpose, we denote the density of the normal distribution corresponding to a regression model (2.1) by $f_{s}(~{}\cdot\mid x_{i},\theta_{s})$ with parameter $\theta_{s}=(\sigma^{2},\vartheta_{s})^{\top}$ and identify the different models by their densities $f_{1},\ldots,f_{r}$ (note that in the situation considered in this sections these only differ in the mean). Using the observations $y_{n}=(y_{11},\ldots,y_{1n_{1}},y_{21},\ldots,y_{kn_{k}})^{\top}$ we calculate in each model the maximum likelihood estimate

[TABLE]

of the parameter $\theta_{s}$, where

[TABLE]

is the log-likelihood in candidate model $f_{s}$ ($s=1,\ldots r$). Note, that we do not assume that the true data generating density is included in the set of candidate models $f_{1},\ldots,f_{r}$. Each estimate $\hat{\theta}_{n,s}$ of the parameter $\theta_{s}$ yields an estimate $\hat{\mu}_{s}=\mu_{s}(\hat{\theta}_{n,s})$ for the quantity of interest, where $\mu_{s}$ is the target parameter in model $s$.

For example, regression models of the type (2.1) are frequently used in dose finding studies (see MacDougall, (2006) or Bretz et al., (2008)). In this case a typical target function $\mu_{s}$ of interest is the “quantile” defined by

[TABLE]

The value defined in (2.4) is well-known as $\text{ED}_{\alpha}$, that is, the effective dose at which $100\times\alpha\%$ of the maximum effect in the design space ${\cal X}=[a,b]$ is achieved.

We now briefly discuss the principle of model selection and averaging to estimate the target parameter $\mu$. For model selection we choose the model $f_{s^{*}}$ from $f_{1},\ldots,f_{s}$, which maximizes Akaike’s information criterion (AIC)

[TABLE]

where $p_{s}$ is the number of parameters in model $f_{s}$ (see Claeskens and Hjort, (2008), Chapter 2). The target parameter is finally estimated by $\hat{\mu}=\mu_{s^{*}}(\hat{\theta}_{n,{s^{*}}})$. Obviously, other model selection schemes, such as the Bayesian or focussed information criterion can be used here as well, but we restrict ourselves to the AIC for the sake of a transparent presentation.

Roughly speaking, model averaging is a weighted average of the individual estimates in the competing models. It might be viewed from a Bayesian (see for example Wassermann, (2000)) or a frequentist point of view (see for example Claeskens and Hjort, (2008)) resulting in different choices of model averaging weights. We will focus here on non-Bayesian methods. More explicitly, assigning nonnegative weights $w_{1},\ldots,w_{r}$ to the candidate models $f_{1},\ldots,f_{r},$ with $\sum_{i=1}^{r}w_{i}=1$, the model averaging estimate for $\mu$ is given by

[TABLE]

Frequently used weights are uniform weights (see, for example Stock and Watson, (2004), Kapetanios et al., (2008)). More elaborate model averaging weights can be chosen depending on the data. For example, Claeskens and Hjort, (2008) define smooth AIC-weights as

[TABLE]

Alternative data dependent weights can be constructed using other information criteria or model selection criteria. There also exists a vast amount of literature on determining optimal data dependent weights such that the resulting mean squared error of the model averaging estimate is minimal (see Hjort and Claeskens, (2003), Hansen, (2007) or Liang et al., (2011) among many others). For the sake of brevity we concentrate on smooth AIC-weights here, but similar observations as presented in this paper can also be made for other data dependent weights.

2.2 The class of competing models matters

In this section we illustrate the influence of the candidate set on the properties of model averaging estimation and estimation after model selection by means of a brief simulation study. For this purpose we consider four regression models of the form (2.1), which are commonly used in dose-response modeling and specified in Table 1 with corresponding parameters.

Here we adapt the setting of Pinheiro et al., (2006) who model the dose-response relationship of an anti-anxiety drug, where the dose of the drug may vary in the interval ${\cal X}=[0,150]$. In particular, we have $k=6$ different dose levels $x_{i}\in\{0,10,25,50,100,150\}$ and patients are allocated to each dose level most equally, where the total sample size is $n\in\{50,100,250\}$. We consider the problem of estimating the $\text{ED}_{0.4}$, as defined in (2.4).

To investigate the particular differences between both estimation methods we choose two different sets of competing models from Table 1. The first set

[TABLE]

contains the log-linear, the Emax and the quadratic model, while the second set

[TABLE]

contains the log-linear, the Emax and the exponential model. The set ${\cal S}_{1}$ serves as a prototype set of “similar” models while the set ${\cal S}_{2}$ contains models of more “different” shape. This is illustrated in Figure 1. In the left panel we show the quadratic model $f_{4}$ (for the parameters specified in Table 1) and the best approximations of this function by a log-linear model ($f_{1}$) and an Emax model ($f_{2}$) with respect to the Kullback-Leibler divergence

[TABLE]

In this case, all models have a very similar shape and we obtain for the ED0.4 the values $32.581$, $32.261$ and $33.810$ for the log-linear ($f_{1}$), Emax ($f_{2}$) and quadratic model ($f_{4}$). Similarly the right panel shows the exponential model ($f_{3}$, solid line) and its corresponding best approximations by the log-linear model ($f_{1}$) and the Emax model ($f_{2}$). Here we observe larger differences between the models in the candidate set and we obtain for the ED0.4 the values $58.116$, $42.857$ and $91.547$ for the models $f_{1}$, $f_{2}$ and $f_{3}$, respectively.

All results presented in this paper are based on $1000$ simulations runs generating in each run $n$ observations of the form

[TABLE]

where the errors $\varepsilon_{ij}^{(l)}$ are independent centered normal distributed random variables with $\sigma^{2}=0.1$ and $\eta_{s}$ is one of the models $\eta_{1},\ldots,\eta_{4}$ (with parameters specified in Table 1). The parameter $\mu=\text{ED}_{0.4}$ is estimated by model averaging with uniform weights, smooth AIC-weights given in (2.7) and estimation after model selection by the AIC criterion.

In Table 2 and 3 we show the simulated mean squared errors of the model averaging estimates with uniform weights (left column), smooth AIC-weights (middle column) and estimation after model selection (right column). Here, different rows correspond to different models. The numbers printed in bold face indicate the estimation method with the smallest mean squared error.

2.2.1 Models of similar shape

We will first discuss the results for the set of similar models in (2.8) (see Table 2). If the data generating model is an element of the set of candidate models, model averaging with uniform weights performs very well. Model averaging with smooth AIC-weights yields an about $10\%$ -$25\%$ larger mean squared error (except for two cases, where it performs better than model averaging with uniform weights). On the other hand the mean squared error of estimation after model selection is substantially larger than that of model averaging, if the sample size is small. This is a consequence of the additional variability associated with data-dependent weights. For example, if the sample size is $n=50$ and the data generating model is given by $f_{1}$, the mean squared errors of the model averaging estimates with uniform and smooth AIC-weights and the estimate after model selection are given by $437.0$, $498.3$ and $759.0$, respectively. The corresponding variances are given by $235.2$, $337.6$ and $599.7$, respectively. For the squared bias the order is exactly the opposite, that is $201.9$, $160.7$, $159.3$, but the differences are not so large. This means that the bias can be reduced by using random weights, because these put more weight on the “correct” model. As a consequence, compared to model averaging with uniform weights the performance of model averaging with smooth AIC-weights and the estimate after model selection improves with increasing sample size. Nevertheless, if the “true” model is an element of the candidate set and the functions in this set have a similar shape, model averaging performs better than estimation after model selection. In particular, model averaging with (fixed) uniform weights yields very reasonable results. These observations coincide with the findings of Schorning et al., (2016) and Buatois et al., (2018) who compared model averaging and model selection in the context of dose finding studies (see also Chen et al., (2018) for similar results for the AIC in the context of ordered probit and nested logit models).

The situation changes if none of the candidate models from the set ${\cal S}_{1}$ is the “true” model. This is illustrated in the lower part of Table 2, where we show results if the exponential model $f_{3}$ is used for generating the data. We observe that model averaging with uniform weights is outperformed by model averaging with smooth AIC-weights. Moreover, the estimate after model selection is even better, if the sample size increases. These observations can be explained by the different shapes of the regression functions, as illustrated in Figure 1. By a suitable choice of parameters the quadratic model can adapt to the shape of the exponential model, whereas the log-linear and the Emax model still have different forms (see right panel of Figure 1 for the best approximations that are possible using the log-linear and the Emax model). Thus, incorporating these models in a model averaging estimate yields a large bias, that can be reduced substantially by data dependent weights or by model selection. For example, if $n=100$ the squared bias of the model averaging estimate with uniform weights is $981.631$, whereas the model averaging estimate with smooth AIC-weights and the estimate after model selection show a squared bias of $328.634$ and $69.465$, respectively.

2.2.2 Models of more different shape

We will now consider the candidate set ${\cal S}_{2}$ in (2.9), which serves as an example of more different models and includes the log-linear, the Emax and the exponential model. The simulated mean squared errors of the three estimates of the ED0.4 are given in Table 3. The upper part of the table corresponds to cases, where data is generated from a model in the candidate set ${\cal S}_{2}$ used for model selection and averaging. In contrast to Section 2.2.1 we observe only one scenario, where model averaging with uniform weights gives the smallest mean squared error (but in this case model averaging with smooth AIC-weights yields very similar results). If the sample size increases model averaging with smooth AIC-weights and estimation after model selection yield a substantially smaller mean squared error. An explanation of this observation consists in the fact that for a candidate set containing models with a rather different shape model averaging with uniform weights produces a large bias. On the other hand model averaging with smooth AIC-weights and estimation after model selection adapt to the data and put more weight on the “true” model, in particular if the sample size is large. As estimation after model selection has a larger variance and the variance is decreasing with increasing sample size, the bias is dominating the mean squared error for large sample sizes and thus estimation in the model selected by the AIC is more efficient for large sample sizes.

Finally, if the data is generated according to the quadratic model $f_{4}\not\in{\cal S}_{2}$ model averaging with uniform weights has the smallest mean squared error if the sample is $n=50$ and $n=100$. In this case estimation in the model selected by the AIC performs much worse (due to its large variance). However, the differences become smaller with increasing sample size. In particular for $n=250$ model averaging with smooth AIC-weights and estimation after model selection show a substantially better performance than model averaging with uniform weights.

The numerical study in Section 2.2.1 and 2.2.2 can be summarized as follows. The results observed in the literature have to be partially relativized. The superiority of model averaging with uniform weights can only be observed for classes of “similar” competing models and a not too large signal to noise ratio. On the other hand if the models in the candidate set are of rather different structure, model averaging with data dependent weights (such as smooth AIC-weights) or estimation after model selection may show a better performance. For these reasons we will investigate optimal/efficient designs for all three estimation methods in the following sections. We will demonstrate that a careful design of experiments can improve the accuracy of these estimates substantially.

3 Asymptotic properties and optimal design

In this section we will derive the asymptotic properties of model averaging estimates with fixed weights in the case where the competing models are not nested. The results can be used for (at least) two purposes. On the one hand they provide some understanding of the empirical findings in Section 2, where we observed, that for increasing sample size the mean squared error of model averaging estimates is dominated by its bias. On the other hand, we will use these results to develop an asymptotic representation of the mean squared error of the model averaging estimate, which can be used in the construction of optimal designs.

3.1 Model averaging for non-nested models

Hjort and Claeskens, (2003) provide an asymptotic distribution of frequentist model averaging estimates making use of local alternatives which require the true data generating process to lie inside a wide parametric model. All candidate models are sub-models of this wide model and the deviations in the parameters are restricted to be of order $n^{-1/2}$. Using this assumption results in convenient approximations for the mean squared error as variance and bias are both of order $O(1/n)$. However, in the discussion of this paper Raftery and Zheng, (2003) pose the question if the framework of local alternatives is realistic. More importantly, frequentist model averaging is also often used for non-nested models (see for example Verrier et al., (2014)). In this section we will develop asymptotic theory for model averaging estimation in non-nested models. In particular, we do not assume that the “true” model is among the candidate models used in the model averaging estimate.

As we will apply our results for the construction of efficient designs for model averaging estimation we use the common notation of this field. To be precise, let $Y$ denote a response variable and let $x$ denote a vector of explanatory variables defined on a given compact design space ${\cal X}$. Suppose that $Y$ has a density $g(y\mid x)$ with respect to a dominating measure. For estimating a quantity of interest, say $\mu$, from the distribution $g$ we use $r$ different parametric candidate models with densities

[TABLE]

where $\theta_{s}$ denotes the parameter in the $s$th model, which varies in a compact parameter space, say $\Theta_{s}\subset\mathbb{R}^{p_{s}}$ $(s=1,...,r)$. Note, that in general we do not assume that the density $g$ is contained in the set of candidate models in (3.1) and that the regression model (2.1) investigated in Section 2 is a special case of this general notation.

We assume that $k$ different experimental conditions, say $x_{1},\ldots,x_{k}$, can be chosen in a design space $\mathcal{X}$ and that at each experimental condition $x_{i}$ one can observe $n_{i}$ responses, say $y_{i1},\ldots,y_{in_{i}}$ (thus the total sample size is $n=\sum_{i=1}^{k}n_{i}$), which are realizations of independent identically distributed random variables $Y_{i1},\ldots,Y_{in_{i}}$ with density $g(\cdot\mid x_{i})$. For example, if $g$ coincides with $f_{s}$ then the density of the random variables $Y_{i1},\ldots,Y_{in_{i}}$ is given by $f_{s}(~{}\cdot\mid x_{i},\theta_{s})$ ($i=1,\ldots,k$). To measure efficiency and to compare different experimental designs we will use asymptotic arguments and consider the case $\lim_{n\to\infty}\tfrac{n_{i}}{n}=\xi_{i}\in(0,1)$ for $i=1,\ldots,k$. As common in optimal design theory we collect this information in the form

[TABLE]

which is called approximate design in the following discussion (see, for example, Kiefer, (1974)). For an approximate design $\xi$ of the form (3.2) and total sample size $n$ a rounding procedure is applied to obtain integers $n_{i}$ taken at each $x_{i}$ ($i=1,\ldots,{k})$ from the not necessarily integer valued quantities $\xi_{i}n$ (see, for example Pukelsheim, (2006), Chapter 12).

The asymptotic properties of the maximum likelihood estimate (calculated under the assumption that $f_{s}$ is the correct density) is derived under certain assumptions of regularity (see the Assumptions (A1)-(A6) in Section 6). In particular, we assume that the functions $f_{s}$ are twice continuously differentiable with respect to $\theta_{s}$ and that several expectations of derivatives of the log-densities exist. For a given approximate design $\xi$ and a candidate density $f_{s}$ we denote by

[TABLE]

the Kullback-Leibler divergence between the models $g$ and $f_{s}$ and assume that

[TABLE]

is unique for each $s\in\{1,\ldots,r\}$. For notational simplicity we will omit the dependency of the minimum on the density $g$, whenever it is clear from the context and denote the minimizer by $\theta_{s}^{*}(\xi)$. We also assume that the matrices

[TABLE]

exist, where expectations are taken with respect to the true distribution ${g(\cdot\mid x_{i})}$.

Under standard assumptions White, (1982) shows the existence of a measurable maximum likelihood estimate $\hat{\theta}_{n,s}$ for all candidate models which is strongly consistent for the (unique) minimizer $\theta_{s}^{*}(\xi)$ in (3.4). Moreover, the estimate is also asymptotically normal distributed, that is

[TABLE]

where we assume the existence of the inverse matrices, $\overset{\mathcal{D}}{\longrightarrow}$ denotes convergence in distribution and we use the notations

[TABLE]

($s,t=1,\ldots r$). The following result gives the asymptotic distribution of model averaging estimates of the form (2.6).

Theorem 3.1**.**

If Assumptions (A1) - (A7) in Section 6 are satisfied, then the model averaging estimate (2.6) satisfies

[TABLE]

where the asymptotic variance is given by

[TABLE]

Theorem 3.1 shows, that the model averaging estimate is biased for the true target parameter $\mu_{\text{true}}$, unless we have $\sum_{s=1}^{r}w_{s}\mu_{s}(\theta_{s}^{*}(\xi))=\mu_{\text{true}}$. Hence we aim to minimize the asymptotic mean squared error of the model averaging estimate. Note, that the bias does not depend on the sample size, while the variance is of order $O(1/n)$.

3.2 Optimal designs for model averaging of non-nested models

Alhorn et al., (2019) determined optimal designs for model averaging minimizing the asymptotic mean squared error of the estimate calculated in a class of nested models under local alternatives and demonstrated that optimal designs lead to substantially more precise model averaging estimates than commonly used designs in dose finding studies. With the results of Section 3.1 we can develop a more general concept of design of experiments for model averaging estimation, which is applicable for non-nested models and in situations, where the “true” model is not contained in the set of candidate models used for model averaging.

To be precise, we consider the criterion

[TABLE]

where $\mu_{\text{true}}$ is the target parameter in the “true” model with density $g$ and $\sigma_{w}^{2}(\theta^{*}(\xi))$ and $\theta_{s}^{*}(\xi)$ are defined in (3.10) and (3.4), respectively. Note that this criterion depends on the “true” distribution via $\mu_{\text{true}}$ and the best approximating parameters $\theta_{s}^{*}(\xi)=\theta_{s,g}^{*}(\xi)$.

For estimating the target parameter $\mu$ via a model averaging estimate of the form (2.6) most precisely a “good” design $\xi$ yields small values of the criterion function $\Phi_{\text{mav}}(\xi,g)$. Therefore, for a given finite set of candidate models $f_{1},\ldots,f_{r}$ and weights $w_{s},s=1,\ldots,r,$ a design $\xi^{*}$ is called locally optimal design for model averaging estimation of the parameter $\mu$, if it minimizes the function $\Phi_{\text{mav}}(\xi,g)$ in (3.11) in the class of all approximate designs on $\mathcal{X}$. Here the term “locally” refers to the seminal paper of Chernoff, (1953) on optimal designs for nonlinear regression models, because the optimality criterion still depends the unkown density $g(y\mid x)$.

A general approach to address this uncertainty problem is a Bayesian approach based on a class of models for the density $g$. To be precise, let ${\cal G}$ denote a finite set of potential densities and let $\pi$ denote a probability distribution on ${\cal G}$, then we call a design Bayesian optimal design for model averaging estimation of the parameter $\mu$ if it minimizes the function

[TABLE]

In general, the set ${\cal G}$ can be constructed independently of the set of candidate models. However, if there is not much prior information available one can construct a class of potential models ${\cal G}$ from the candidate set as follows. We denote the candidate set of models in (3.1) by $\mathcal{S}$. Each of these models depends on a unknown parameter $\theta_{s}$ and we denote by ${\cal F}_{f_{s}}\subset\Theta_{s}$ a set of possible parameter values for the model $f_{s}$. Now let $\pi_{2}$ denote a prior distribution on ${\cal S}$ and for each $f_{s}\in{\cal S}$ let $\pi_{1}(\cdot\mid f_{s})$ denote a prior distribution on ${\cal F}_{f_{s}}$. Finally, we define $\mathcal{G}=\{(g,\theta):g\in\mathcal{S},\theta\in\mathcal{F}_{g}\}$ and a prior

[TABLE]

then the criterion (3.12) can be rewritten as

[TABLE]

In the finite sample study of the following section the set ${\cal S}$ and the set $\mathcal{F}_{g}$ (for any $g\in{\cal S}$) are finite, which results in a finite set ${\cal G}$.

Locally and Bayesian optimal designs for model averaging estimation have to be calculated numerically in all cases of practical interest. We will state now a necessary condition for the optimality of a given design with respect to the criterion $\phi_{\text{mav}}^{\pi}$. Note, that this criterion is not convex and therefore a sufficient condition cannot be derived. In the following discussion we denote by $A_{s}^{*}=A_{s}(\theta_{s,g}^{*}(\xi^{*}),\xi^{*})$ and $B_{st}^{*}=B_{st}(\theta_{s,g}^{*}(\xi^{*}),\theta_{t,g}^{*}(\xi^{*}),\xi^{*})$ the matrices defined in (3.5) and (3.6), respectively, evaluated in $\xi^{*}$ and $\theta_{s,g}(\xi^{*})$.

Theorem 3.2**.**

If a design $\xi^{*}$ is Bayesian optimal for model averaging estimation of the parameter $\mu$ with respect to the prior $\pi$, then

[TABLE]

*holds for all $x\in\mathcal{X}$, where the derivatives $\theta_{s}^{\prime}(\xi^{*},x)$ and $\sigma_{g}^{\prime}(\xi^{*},x)$ are given by *

[TABLE]

*respectively. Here the matrices $h_{st,g}^{\prime}(\xi^{*},x)$ and $h_{s,g}^{\prime}(\xi^{*},x)$ are given by *

[TABLE]

where the matrix

[TABLE]

contains the third derivatives of the log-likelihood with respect to the parameters $\theta_{s}=(\theta_{s,1},\ldots,\theta_{s,p_{s}})^{\top}$. Moreover, there is equality in (3.2) for all support points of the optimal design.

Example 3.1**.**

We illustrate the application of Theorem 3.2 for regression models of the from (2.1) with centred normal distributed errors. As regression functions we use the log-linear and the Emax model and their parameter specifications given in Table 1. Then, the locally optimal design for estimation of the $\text{ED}_{0.4}$ in the log-linear model $f_{1}$ and in the Emax model $f_{2}$ are given by

[TABLE]

and $\left\{0,18.75,150;~{}0.25,0.5,0.25\right\},$ respectively [see Dette et al., (2010)]. For sample size $n=100$ we determine a Bayesian optimal design for model averaging estimation of the $\text{ED}_{0.4}$ (with uniform weights) with respect to the criterion (3.12). The set of possible models is given by $\mathcal{G}=\{f_{1},f_{2}\}$ with parameters specified in Table 1, and we choose a uniform prior on this set. The optimal design has been calculated numerically using the COBYLA algorithm (see Powell, (1994)) and is given by

[TABLE]

The necessary condition of Theorem 3.2 is satisfied as illustrated in Figure 2). Note that the design $\xi_{12}^{*}$ can be considered as a compromise between the locally optimal designs for the individual models and that $\xi_{12}^{*}$ would not be optimal if the inequality was not satisfied. **

We conclude noting that the optimality criteria proposed in this section have been derived for model averaging estimates with fixed weights. The asymptotic theory presented here cannot be easily adapted to estimates using data-dependent (random) weights (as considered in Section 2), because it is difficult to get an explicit expression for the asymptotic distribution, which is not normal in general. Nevertheless, we will demonstrate in the following section that designs minimizing the mean squared error of model averaging estimates with fixed weights will also yield a substantial improvement in model averaging estimation with smooth AIC-weights and in estimation after model selection.

4 Bayesian optimal designs for model averaging

We will demonstrate by means of a simulation study that the performance of all considered estimates can be improved substantially by the choice of an appropriate design. For this purpose we consider the same situation as in Section 2, that is regression models of the from (2.1) with centred normal distributed errors. We also consider the two different candidate sets ${\cal S}_{1}$ and ${\cal S}_{2}$ defined in (2.8) (log-linear, Emax and quadratic model) and (2.9) (log-linear, Emax and exponential model), respectively.

Using the criterion introduced in Section 3 we now determine a Bayesian optimal design for model averaging estimation of the $\text{ED}_{0.4}$ with uniform weights from $n=100$ observations. We require a prior distribution for the unknown density $g$, and we use a distribution of the form (3.13) for this purpose. To be precise, let $f_{s}(y\mid x,\theta_{s})$ denote the density of a normal distribution with mean $\eta_{s}(x,\vartheta_{s})$ and variance $\sigma^{2}_{s}=0.1$ ($s=1,\ldots,r$), where the mean functions are given in Table 1. As the criterion (3.14) does not depend on the intercept $\vartheta_{s1}$, these are not varied and taken from Table 1. For each of the other parameters we use three different values: the values specified in Table 1 and a $10\%$ larger and smaller value of this parameter.

[TABLE]

4.1 Models of similar shape

We will first consider the candidate set ${\cal S}_{1}=\{f_{1},f_{2},f_{4}\}$ consisting of the log-linear, the Emax and the quadratic model. For the definition of the prior distribution (3.13) in the criterion (3.14) we consider a uniform distribution $\pi_{2}$ on the set ${\cal S}_{1}$ and a uniform prior $\pi_{1}(\cdot\mid f_{s})$ on each set ${\cal F}_{f_{s}}$ in (4.1) ($s=1,2,4$). The Bayesian optimal design for model averaging estimation of the $\text{ED}_{0.4}$ minimizing the criterion (3.14) is given by

[TABLE]

We will compare this design with the design

[TABLE]

proposed in Pinheiro et al., (2006) for a a similar setting (this design has also been used in Section 2) and the locally optimal design for the estimation of the $\text{ED}_{0.4}$ in the log-linear model given by (3.25). Results for the locally optimal designs for estimation of the $\text{ED}_{0.4}$ in the Emax and exponential model are similar and omitted for the sake of brevity. We use the same setup as in Section 2. Only results for the sample size $n=100$ are reported and results for other sample sizes are available from the authors.

The corresponding results are given in Table 4, where we use the models $f_{1},f_{2},f_{3}$ and $f_{4}$ from Table 1 to generate the data (note that the model $f_{3}$ is not in the candidate set used for model averaging and model selection). The different columns represent the different estimation methods (left column: model averaging with uniform weights; middle column: smooth AIC-weights, right column: model selection). The numbers printed in boldface indicate the minimal mean squared error for each estimation method obtained from the different experimental designs. First, we consider the situation, where the data generating model is contained in the set of candidate models ${\cal S}_{1}=\{f_{1},f_{2},f_{4}\}$ corresponding to the upper part of the table. We observe that in this case model averaging yields better results than estimation after model selection and this superiority is independent of the design under consideration. Compared to the designs $\xi_{1}$ and $\xi_{2}$ the Bayesian optimal design $\xi_{{\cal S}_{1}}^{*}$ for model averaging with uniform weights improves the efficiency of all estimation techniques. For example, when data is generated using the log-linear model $f_{1}$ the mean squared error of the model averaging estimate with uniform weights is reduced by $20.5\%$ and $4.2\%$, when the optimal design is used instead of the designs $\xi_{1}$ or $\xi_{2}$, respectively. This improvement is remarkable as the design $\xi_{2}$ is locally optimal for estimating the ED0.4 in the model $f_{1}$ and data is generated from this model. In other cases the improvement is even more visible. For example, if data is generated by the model $f_{2}$ the improvement in model averaging estimation with uniform weights is $25.1\%$ and $71.7\%$ compared to the designs $\xi_{1}$ and $\xi_{2}$ defined in (4.3) and (3.25). Moreover, although the designs are constructed for model averaging with uniform weights they also yield substantially more accurate model averaging estimates with smooth AIC-weights and a more precise estimate after model selection. For example, if the data is generated from model $f_{1}$ the mean squared error is reduced by $24.2\%$ and by $10.2\%$ for estimation with smooth AIC-weights and by $39.1\%$ and $49.1\%$ for estimation after model selection, respectively. Similar results can be observed for the models $f_{2}$ and $f_{4}$.

Next, we consider the case where the data is generated from the exponential model $f_{3}$, which is not contained in the candidate set ${\cal S}_{1}$. The efficiency of all three estimates improves substantially by the use of the Bayesian optimal design $\xi_{{\cal S}_{1}}^{*}$. Interestingly, the improvement is less pronounced for model averaging with uniform weights ($4.6\%$ and $66.8\%$ compared to the designs $\xi_{1}$ and $\xi_{2}$ in in (4.3) and (3.25), respectively) than for smooth AIC-weights ($10.5\%$ and $77.5\%$) and estimation after model selection ($16.9\%$ and $85.2\%$).

Summarizing, our numerical results show that the Bayesian optimal design for model averaging estimation of the $\text{ED}_{0.4}$ yields a substantial improvement of the mean squared error of the model averaging estimate with uniform weights ($4.2\%$-$71.7\%$), smooth AIC-weights ($10.2\%$-$77.5\%$) and the estimate after model selection ($16.9\%$-$85.4\%$) for all four models under consideration.

4.2 Models of different shape

We will now consider the second candidate set ${\cal S}_{2}$ consisting of the log-linear ($f_{1}$) the Emax ($f_{2}$) and the exponential model ($f_{3}$). For the definition of the prior distribution (3.13) in the criterion (3.14) we use a uniform distribution $\pi_{2}$ on the set ${\cal S}_{2}$ and a uniform prior $\pi_{1}(\cdot\mid f_{s})$ on each set ${\cal F}_{f_{s}}$ ($s=1,2,3$) in (4.1). For this choice the Bayesian optimal design for model averaging estimation of the $\text{ED}_{0.4}$ is given by

[TABLE]

and has (in comparison to the design $\xi_{{\cal S}_{1}}^{*}$ in Section 4.1) five instead of four support points.

The simulated mean squared errors of the three estimates under different designs are given in Table 5. We observe again that compared to the designs $\xi_{1}$ and $\xi_{2}$ in (4.3) and (3.25) the Bayesian optimal design $\xi^{*}_{{\cal S}_{2}}$ improves most estimation techniques substantially. However, if model averaging with uniform weights is used and data is generated by model $f_{2}$ or $f_{3}$, the mean squared error of the model averaging estimate from the optimal design is $5.4\%$ and $4.5\%$ larger than the mean squared error obtained by the design $\xi_{1}$, respectively. For model averaging with smooth AIC-weights and data being generated from model $f_{2}$ this difference is $5.9\%$. Overall, the reported results demonstrate a substantial improvement in efficiency by usage of the Bayesian optimal design independently of the estimation method. If the Bayesian optimal design is used, estimation after model selection yields the smallest mean squared error if the data is generated from a model of the candidate set ${\cal S}_{2}$. On the other hand, if data is generated from model $f_{4}\not\in{\cal S}_{2}$ model averaging with equal weights shows the best performance.

Summarizing, our numerical results show that compared to the designs $\xi_{1}$ and $\xi_{2}$ the design $\xi^{*}_{{\cal S}_{2}}$ reduces the mean squared error of model averaging estimates with uniform weights up to $69.4\%$. Furthermore, for smooth AIC-weights and estimation after model selection the reduction can be even larger and is up to $70.5\%$ and $85.3\%$, respectively. These improvements hold also for the quadratic model $f_{4}$, which is not contained in the candidate set ${\cal S}_{2}$ used in the definition of the optimality criterion.

5 Conclusions

In this paper we derived the asymptotic distribution of the frequentist model averaging estimate with fixed weights from a class of not necessarily nested models. We neither assume that this class contains the “true” model. We use these results to determine Bayesian optimal designs for model averaging, which can improve the estimation accuracy of the estimate substantially. Although these designs are constructed for model averaging with fixed weights, they also yield a substantial improvement of accuracy for model averaging with data dependent weights and for estimation after model selection.

We also demonstrate that the superiority of model averaging against estimation after model selection depends sensitively on the class of competing models, which is used in the model averaging procedure. If the competing models are similar (which means that a given model from the class can be well approximated by all other models), then model averaging should be preferred. Otherwise, we observe advantages for estimation after model selection, in particular, if the signal to noise ratio is small.

Although, the new designs show a very good performance for estimation after model selection and for model averaging with data dependent weights, it is of interest to develop optimal designs, which address the specific issues of data dependent weights in the estimates. This is a very challenging problem for future research as there is no simple expression of the asymptotic mean squared error of these estimates. A first approach to solve this problem is an adaptive one and a further interesting and very challenging question of future research is to improve the accuracy of adaptive designs.

6 Technical assumptions and proofs

6.1 Assumptions

Following White, (1982) we assume:

• (A1)

The random variables $Y_{ij},i=1,\ldots,k,j=1,\ldots,n_{i}$ are independent. Furthermore, $Y_{i1},\ldots,Y_{in_{i}}$ have a common distribution function with a measurable density $g(\cdot\mid x_{i})$ with respect to a dominating measure $\nu$.

• (A2)

The distribution function of each candidate model $s\in\{1,\ldots,r\}$ has a measurable density $f_{s}(\cdot\mid x,\theta_{s})$ with respect to $\nu$ (for all $\theta_{s}\in\Theta_{s}$) that is continuous in $\theta_{s}$.

• (A3)

For all $x\in\mathcal{X}$ the expectation $\mathbb{E}(\log(g(Y\mid x)))$ exists (where expectation is taken with respect to $g(\;\cdot\mid x)$) and for each candidate model the function $y\mapsto|\log f_{s}(y\mid x,\theta_{s})|$ is dominated by a function that is integrable with respect to $g(\;\cdot\mid x)$ and does not depend on $\theta_{s}$. Furthermore the Kullback-Leibler divergence (3.3) has a unique minimum $\theta_{s,g}^{*}(\xi)$ defined in (3.4) and $\theta_{s,g}^{*}(\xi)$ is an interior point of $\Theta_{s}$.

• (A4)

For all $x\in\mathcal{X}$ the function $y\mapsto\tfrac{\partial\log f_{s}(y\mid x,\theta_{s})}{\partial\theta_{s}}$ is a measurable function for all $\theta_{s}\in\Theta_{s}$ and continuously differentiable with respect to $\theta_{s}$ for all $y\in\mathbb{R}$.

• (A5)

The entries of the (matrix valued) functions $\tfrac{\partial^{2}\log f_{s}(y\mid x,\theta_{s})}{\partial\theta_{s}\partial\theta_{s}^{\top}}$, \tfrac{\partial\log f_{s}(y\mid x,\theta_{s})}{\partial\theta_{s}}\big{(}\tfrac{\partial\log f_{t}(y\mid x,\theta_{t})}{\partial\theta_{t}}\big{)}^{\top} are dominated by integrable functions with respect to $g(\;\cdot\mid x)$ for all $x\in\mathcal{X}$ and $\theta_{s}\in\Theta_{s}$.

• (A6)

The matrices $B_{ss}(\theta_{s}^{*}(\xi),\theta_{s}^{*}(\xi),\xi)$ and $A_{s}(\theta_{s}^{*}(\xi),\xi)$ in (3.5) and (3.6) are nonsingular.

• (A7)

The functions $\theta_{s}\mapsto\mu_{s}(\theta_{s})$ are once continuously differentiable.

6.2 Proof of Theorem 3.1.

By equation (A.2) in White, (1982) we have

[TABLE]

where $\overset{p}{\longrightarrow}$ denotes convergence in probability (note that the matrix $A_{s}(\theta_{s}^{*})=A_{s}(\theta_{s}^{*},\xi)$ is nonsingular by assumption). An application of the multivariate central limit theorem now leads to

[TABLE]

where $B_{{st}}=B_{{st}}(\theta_{s}^{*}(\xi),\theta_{t}^{*}(\xi),\xi)$ is defined in (3.6). Combining (A.1) and (A.8) we obtain the weak convergence of the vector $\hat{\theta}_{n}=(\hat{\theta}_{n,1}^{\top},\ldots,\hat{\theta}_{n,r}^{\top})^{\top}$, that is $\sqrt{n}(\hat{\theta}_{n}-\theta^{*}(\xi))\overset{\mathcal{D}}{\longrightarrow}\mathcal{N}(0,\Sigma),$ where $\Sigma=(\Sigma_{st})_{s,t=1,\ldots,r}$ is a block matrix with entries $\Sigma_{{st}}=A_{s}^{-1}(\theta_{s}^{*}(\xi))B_{st}(\theta_{s}^{*}(\xi),\theta_{t}^{*}(\xi))A_{t}^{-1}(\theta_{t}^{*}(\xi))$ ($s,t=1,\ldots,r$) and the vector $\theta_{s}^{*}(\xi)$ is given by $\theta^{*}(\xi)=(\theta_{1}^{*}(\xi)^{\top},\ldots,\theta_{r}^{*}(\xi)^{\top})^{\top}$.

Next, we define for the parameter vector $\theta^{\top}=(\theta_{1}^{\top},...,\theta_{r}^{\top})\in\mathbb{R}^{\sum_{s=1}^{r}p_{s}}$ the projection $\pi_{s}$ by $\pi_{s}\theta:=\theta_{s}$ and the vector \tilde{\mu}(\theta)=\big{(}\mu_{1}(\pi_{1}\theta),\ldots,\mu_{r}(\pi_{r}\theta)\big{)}^{T} with derivative

[TABLE]

An application of the Delta method shows that \sqrt{n}(\tilde{\mu}(\hat{\theta}_{n})-\tilde{\mu}(\theta^{*}(\xi)))\overset{\mathcal{D}}{\longrightarrow}\mathcal{N}\big{(}0,\mu_{\theta^{*}(\xi)}^{\prime}\Sigma(\mu_{\theta^{*}(\xi)}^{\prime})^{\top}\big{)}. The assertion finally follows from the continuous mapping theorem observing the representation $\hat{\mu}_{\text{mav}}=\left(w_{1},\ldots,w_{r}\right)\tilde{\mu}(\hat{\theta}_{n}).$

6.3 Proof of Theorem 3.2.

Throughout this proof we assume that integration and differentiation are interchangeable. Following the arguments in Pukelsheim, (2006), Chapter 11, a Bayesian optimal design $\xi^{*}$ for model averaging estimation of the parameter $\mu$ satisfies the inequality

[TABLE]

for all $x\in\mathcal{X}$, where $D\Phi_{\text{mav}}(\xi^{*},g)(\xi_{x}-\xi^{*})$ denotes the directional derivative of the function $\Phi_{\text{mav}}$ evaluated in the optimal design $\xi^{*}$ in direction $\xi_{x}-\xi^{*}$ and $\xi_{x}$ denotes the Dirac measure at the point $x\in{\cal X}$.

To calculate the derivative we start with the derivative of the parameter $\theta_{s,g}^{*}(\xi)$ defined in (3.4) and define $\theta_{s,g}(\alpha):=\theta_{s,g}^{*}(\xi_{\alpha})$ for $\xi_{\alpha}=\alpha\xi_{x}+(1-\alpha)\xi^{*}$. Note that $\theta_{s,g}(\alpha)$ is the solution of the equation

[TABLE]

and that the derivatives of the left hand side are given by

[TABLE]

By the implicit function theorem we get \tfrac{\partial\theta_{s,g}(\alpha)}{\partial\alpha}=-\big{(}\tfrac{\partial F_{s,g}}{\partial\theta_{s}}\big{)}^{-1}\tfrac{\partial F_{s,g}}{\partial\alpha} and hence

[TABLE]

where $\theta_{s}^{\prime}(\xi^{*},x)$ is defined in (3.2). Consider now the directional derivative of the matrix $B_{st}$ defined in (3.6). An application of chain and product rule gives

[TABLE]

where $h_{st,g}^{\prime}(\xi^{*},x)$ is defined in (3.2). In a similar way the derivative of the matrix $A_{s}$ defined in (3.5) can be determined. First, using the chain rule, we observe that with $\theta_{s,g}^{\prime}(\xi^{*},x)=(\theta_{s,g,1}^{\prime}(\xi^{*},x),\cdots,\theta_{s,g,p_{s}}^{\prime}(\xi^{*},x))^{\top}$

[TABLE]

where $D_{s}$ is defined in Theorem 3.2. We now observe, that \tfrac{\partial A_{s}(\theta_{s,g}^{*}(\xi_{\alpha}),\xi_{\alpha})}{\partial\alpha}\big{|}_{\alpha=0}=h_{s,g}^{\prime}(\xi^{*},x), where $h_{s,g}^{\prime}(\xi^{*},x)$ is defined in (3.2). Noting, that

[TABLE]

equation (3.17) results by an application of the product rule and combination of the derivatives given above. Finally, we have

[TABLE]

and (3.2) follows.

The proof that there is equality in (3.2) for all support points of the optimal design $\xi^{*}$ follows by a standard argument and the details are omitted for the sake of brevity.

Acknowledgements This work has also been supported in part by the Collaborative Research Center “Statistical modeling of nonlinear dynamic processes” (SFB 823, Teilprojekt C2, T1) of the German Research Foundation (DFG).