Consider the standard linear regression model \(y=f(x)'\theta + \varepsilon\), \(\theta \in \mathbb{R}^m\). For different observations, the errors \(\varepsilon\) are assumed to be independent and homoscedastic (their variance is always \(\sigma^2\)). The space \(\mathcal{X}\) of candidate design points is assumed to be rich enough to guarantee that \(\{f(x): x \in \mathcal{X}\}\) spans \(\mathbb{R}^m\).

For a finite \(\mathcal{X}\), the Elfving theorem suggests that the problem of c-optimality can be reduced to linear programming (https://doi.org/10.1016/j.csda.2008.06.023). Namely, consider the following linear programming problem \[ \left.\begin{array}{rl} \hbox{max} & \alpha\\ \hbox{s.t.} & \begin{bmatrix} F & -h \\ 1'_{2n} & 0 \end{bmatrix} \begin{bmatrix} u \\ \alpha \end{bmatrix} = \begin{bmatrix} 0_m\\ 1 \end{bmatrix}\\ & \begin{bmatrix} u \\ \alpha \end{bmatrix}\geq 0, \end{array}\right. \] where \(F=(f(x_1),\ldots,f(x_n),-f(x_1),\ldots,-f(x_n))\).

Then \(w\) is a c-optimal design and \(\alpha^{-2}\) is the c-optimal variance iff \(w_i=u_i+u_{i+n}\) for the LP solution \((u',\alpha)'\).

This method is implemented in the procedure od_REX of OptimalDesign.

Example: trigonometric model

Compute the approximate c-optimal designs for estimation of the parameters \(\theta_1\) and \(\theta_2\) in the trigonometric model \[ y(x)=\cos(x)\theta_1+\sin(x)\theta_2 + \varepsilon(x), \] where \(\mathcal{X}=[0,\pi]\) is discretized into \(101\) equidistant points.

Ft <- Fx_cube(~I(cos(x1))+I(sin(x1))-1, 0, pi, 101, echo = FALSE)
res <- od_REX(Ft, crit = "c", h = c(1, 0), echo = FALSE)
## Linear programming reformulation will be used for c-optimality.
## [1] Method od_c_LP finished after 0.02 seconds at 2021-05-16 14:04:18 .
## [1] c-criterion value: 1
rbind(res$supp, res$w.supp)
##      [,1]
## [1,]    1
## [2,]    1
res2 <- od_REX(Ft, crit = "c", h = c(0, 1), echo = FALSE)
## Linear programming reformulation will be used for c-optimality.
## [1] Method od_c_LP finished after 0.01 seconds at 2021-05-16 14:04:18 .
## [1] c-criterion value: 1
rbind(res2$supp, res2$w.supp)
##      [,1]
## [1,]   51
## [2,]    1
par(mfrow=c(1,2))
Elf2(Ft, vec=c(0,1))
Elf2(Ft, vec=c(1,0))

Optimization problems leading to c-optimality

Besides designs that minimize the variance of the single parameters, c-optimality is useful in many more practical problems.

Example: quadratic regression

Consider the model \[ y(x)=\theta_0+\theta_1x+\theta_2x^2+\varepsilon(x), \] where \(x\in\mathcal{X}=[0,1]\).

Fq <- Fx_cube(~x1+I(x1^2), 0, 1, 2000, echo=FALSE)

Let’s compute the optimal designs that minimize the estimate of the following quantities:

  1. The mean value of the response in the point 1/2

\[ y\left(\frac{1}{2}\right)=\theta_0+\frac{\theta_1}{2}+\frac{\theta_2}{4} \]

That is, we wish to minimize the variance of the estimate of the linear combination \(h^T\theta\), where \(h^T=(1,\frac{1}{2},\frac{1}{4})\).

c1 <- od_REX(Fq, crit="c", h=c(1, 1/2, 1/4), echo=FALSE)
## Linear programming reformulation will be used for c-optimality.
## [1] Method od_c_LP finished after 0.01 seconds at 2021-05-16 14:04:18 .
## [1] c-criterion value: 1.31249868596337
od_plot(Fq, c1$w.best, X = Fq[,2:3], asp = 1, echo = FALSE, dd.color="red", dd.size=1, dd.pch=16)
## Warning in dirder(Fx, w, crit = crit, h = h, echo = FALSE): The information
## matrix is badly conditioned and there may be problems computing the directional
## derivative.

## $call
## od_plot(Fx = Fq, w = c1$w.best, X = Fq[, 2:3], dd.color = "red", 
##     dd.size = 1, dd.pch = 16, asp = 1, echo = FALSE)
## 
## $Pools
## [1] "Pool return suppressed."
Elf3(Fq, vec=c(1, 1/2, 1/4))
  1. The area under the response curve (AUC)

\[ A=\int_0^1(\theta_0+\theta_1x+\theta_2x^2)dx=\theta_0+\frac{\theta_1}{2}+\frac{\theta_2}{3}, \] hence, \(h^T=(1,\frac{1}{2},\frac{1}{3})\).

c2 <- od_REX(Fq, crit="c", h=c(1, 1/2, 1/3), echo=FALSE)
## Linear programming reformulation will be used for c-optimality.
## [1] Method od_c_LP finished after 0 seconds at 2021-05-16 14:04:20 .
## [1] c-criterion value: 1.36111111111111
od_plot(Fq, c2$w.best, X = Fq[,2:3], asp = 1, echo = FALSE, dd.color="red", dd.size=1, dd.pch=16)

## $call
## od_plot(Fx = Fq, w = c2$w.best, X = Fq[, 2:3], dd.color = "red", 
##     dd.size = 1, dd.pch = 16, asp = 1, echo = FALSE)
## 
## $Pools
## [1] "Pool return suppressed."
Elf3(Fq, vec=c(1, 1/2, 1/3))
  1. The slope of the tangent to the curve \(y(x)=\theta_0+\theta_1x+\theta_2x^2\) in the point \(x=1\)

\[ \frac{d}{dx}(\theta_0+\theta_1x+\theta_2x^2)|_{x=1}=(\theta_1+2\theta_2x)|_{x=1}=\theta_1+2\theta_2 \] with \(h^T=(0,1,2)\).

c3 <- od_REX(Fq, crit="c", h=c(0, 1, 2), echo=FALSE)
## Linear programming reformulation will be used for c-optimality.
## [1] Method od_c_LP finished after 0.02 seconds at 2021-05-16 14:04:22 .
## [1] c-criterion value: 0.0781249608984129
od_plot(Fq, c3$w.best, X = Fq[,2:3], asp = 1, echo = FALSE, dd.color="red", dd.size=1, dd.pch=16)

## $call
## od_plot(Fx = Fq, w = c3$w.best, X = Fq[, 2:3], dd.color = "red", 
##     dd.size = 1, dd.pch = 16, asp = 1, echo = FALSE)
## 
## $Pools
## [1] "Pool return suppressed."
Elf3(Fq, vec=c(0, 1, 2))
  1. The optimal prediction of the mean value of the response in the point \(x=3/2\), that is, outside the design space

\[ y\left(\frac{3}{2}\right)=\theta_0+\frac{3\theta_1}{2}+\frac{9\theta_2}{4} \] and \(h^T=(1,\frac{3}{2},\frac{9}{4})\).

c4 <- od_REX(Fq, crit="c", h=c(1, 3/2, 9/4), echo=FALSE)
## Linear programming reformulation will be used for c-optimality.
## [1] Method od_c_LP finished after 0.02 seconds at 2021-05-16 14:04:23 .
## [1] c-criterion value: 0.169642784366023
od_plot(Fq, c4$w.best, X = Fq[,2:3], asp = 1, echo = FALSE, dd.color="red", dd.size=1, dd.pch=16)

## $call
## od_plot(Fx = Fq, w = c4$w.best, X = Fq[, 2:3], dd.color = "red", 
##     dd.size = 1, dd.pch = 16, asp = 1, echo = FALSE)
## 
## $Pools
## [1] "Pool return suppressed."
Elf3(Fq, vec=c(1, 3/2, 9/4))