These questions provide an opportunity to review concepts from the course and to practise the application of Design of Experiments methods in R.
The function referred to in Question 1 below is available in the apts.doe package, which can be installed from GitHub. Make sure you have the latest version of this package before attempting the assessment.
library(devtools)
install_github("statsdavew/apts.doe", quiet = T)
library(apts.doe)1. Fractional factorial designs This question is
based on a real experiment performed by Tribologists at
Southampton.
Engineers in the Materials Science group at the University of
Southampton wish to run an experiment to investigate the impact of
contamination on material corosion in a “pin-on-disc” tribometer (https://en.wikipedia.org/wiki/Tribometer).
There are six factors that can be varied, each of which can take two values:
The experimenters can afford \(n=16\) runs. Prior knowledge from previous experiments suggests the following interactions may be important:
Design an appropriate fractional factorial experiment.
The function material.response simulates data from this experiment. As input, it takes a vector of six values, each either -1 or +1, representing a single run of the experiment. Use this function to generate data for your chosen design, and perform an approriate analysis. What conclusions can you draw?
Note: the Effects function will not work with lm objects which are over-parameterised (ie which contain missing (NA) coefficients). You will therefore need to make sure your lm object is not over-parameterised to draw main and interaction effect plots. As the design is orthogonal (in the non-aliased effects), the easiest way to do this is to refit a linear model in just the effects of interest, eg to plot the interaction between A and B:
materials.lmAB <- lm(y ~ A*B, data = materials.frame)Of course, this doesn’t break the aliasing in the design! So in this case, AB is probably aliased with another interaction (depending on the design you have chosen).
2. Bayesian \(D\)-optimal
design using optim
Adapt the code from the notes to find locally \(D\)-optimal designs with \(n=2\) runs for the Poisson regression model
with linear predictor \(\eta(x) = \beta_0 +
\beta_1 x\) with the log link (log-linear regression).
Fix \(\beta_0=0\) and find designs for \(\beta_1 = 1, \ldots, 6\), restricting to \(-1\le x\le 1\). For each design, plot the design points along with the expected response. There is a pattern in the optimal designs, can you identify it?
The acebayes package has a function aceglm that can find Bayesian optimal designs for generalised linear models. Read the help file, and the use the function to find a Bayesian \(D\)-optimal design when we assume \(\beta_0 = 0\) and a uniform prior distribution for \(\beta_1\) with support \([1,6]\).
Important: when using aceglm set \(N2 = 0\) (otherwise it will result in a nasty crash!).
3. Bayesian \(D\)-optimal
design using acebayes
In the notes, we examine a compartmental model. A simplified version of
this model is given by \[
y_i \sim N\left(\mu(x_i\,;\,\boldsymbol{\theta}), \sigma^2\right)\,,
\] with \[
\mu(x_i\,;\,\boldsymbol{\theta}) = \theta_3\left[\exp\left(-\theta_1
x\right) - \exp\left(-\theta_2 x\right)\right]\,.
\] Use the acenlm function in acebayes to find
a Bayesian \(D\)-optimal design with
\(n=18\) runs, assuming \(x \in [0, 24]\). Assume uniform prior
distributions for \(\theta_1,\theta_2\): \[
\theta_1 \sim U(0.01884, 0.09884);\qquad \theta_2 \sim U(0.298, 8.298)
\] with \(\theta_3 = 21.8\); see
Gotwalt, Jones, and Steinberg (2009) and
Overstall, Woods, and Adamou (2017). Plot
the optimal design and realisations from the model.