In order to see the difference of treatment effects, we must control the nuisance factors if they are expected to be a major source of variation. If the values of the nuisance factors can be measured in advance of the experiment or controlled during the experiment, then they can be taken into account by using blocking factors. There are situations when some nuisance factors or variables may affect the response and they are also measured. For example, we would like to compare the effects of four diets on the weights of month-old piglets. The response variable is the weight three months after starting the diets, which is likely to be affected by the weight at the beginning of the experiment. We would use the starting weight as a covariate. Therefore the employment of covariates is another way to deal with nuance factors or variables. The analysis of covariance of the study to compare treatment effects in the presence of covariances. Note in terms of design in the example, it is a completely randomized design. What is new here we have additional measurements on the experimental units (piglets’ starting weights).
Covariates are often continuous variables. Their effects on the response variable are assumed to follow some parametric form such as linear.
For a completely randomized design with one covariate, the simplest model is
\[\begin{equation} Y_{i j}=\mu+\tau_i+\beta x_{i j}+\epsilon_{i j} \tag{1} \end{equation}\] where \(Y_{i j}\) is the response variable of the \(j\)th unit from the \(i\)th treatment group; \(\tau_i\) is the effect of the \(i\)th treatment, \(x_{i j}\) is the value of the covariate of the \(j\)th unit in the \(i\)th treatment group, \(\epsilon_{i j}\) are i.i.d. \(N\left(0, \sigma^2\right)\).
Here we assume that the effect of the covariate on the response is linear. It can be extended to more complex cases, e.g., polynomial effects.
Estimate the parameters by minimizing \[ \sum_{i=1}^v \sum_{t=1}^{r_i} e_{i t}^2=\sum_{i=1}^v \sum_{t=1}^{r_i}\left(y_{i t}-\mu-\tau_i-\beta\left(x_{i t}-\bar{x}_{. .}\right)\right)^2 . \] Taking derivatives with respect to \(\mu\), \(\tau_i\) and \(\beta\), we get
\[\begin{eqnarray*} && y_{. .}=n \hat{\mu}+\sum_{i=1}^v r_i \hat{\tau}_i \\ && y_{i .}=r_i\left(\hat{\mu}+\hat{\tau}_i\right)+\hat{\beta} \sum_{t=1}^{r_i}\left(x_{i t}-\bar{x}_{. .}\right), \quad i=1, \ldots, v \\ &&\sum_{i=1}^v \sum_{t=1}^{r_i} y_{i t}\left(x_{i t}-\bar{x}_{. .}\right)= \sum_{i=1}^v \sum_{t=1}^{r_i}\left(\hat{\mu}+\hat{\tau}_i\right)\left(x_{i t}-\bar{x}_{. .}\right) +\sum_{i=1}^v \sum_{t=1}^{r_i} \hat{\beta}\left(x_{i t}-\bar{x}_{. .}\right)^2 \end{eqnarray*}\]
subject to \(\sum_{i=1}^v \tau_i=0\).
Solving for these equations to get the least squares estimators for the parameters. Note this estimators are linear functions of \(Y_{it}\) and therefore have normal distributions.
For a completely randomized design and analysis of covariance model (1), a one-way analysis of covariance is used to test the null hypothesis \(H_0:\left\{\tau_1=\tau_2=\cdots=\tau_v\right\}\) against the alternative hypothesis \(H_A\) that at least two of the \(\tau_i\) ’s differ. As with other models, the test is based on the comparison of error sums of squares under the full and reduced models. If the null hypothesis is true with common treatment effect \(\tau_i=0\), then the reduced model is \[ Y_{i t}=\mu+\beta\left(x_{i t}-\bar{x}_{. .}\right)+\epsilon_{i t}. \]
Hence \[ SST=SSE_0-SSE \]
It can be shown that under \(H_0\), \(SST/\sigma^2\) has a \(\chi^2\) distribution with \((v-1)\) degrees of freedom, and \((v-1)\) is the number of parameters the reduced model reduced by.
If the null hypothesis is true, then \[ \operatorname{MST}/\text{MSE} \sim F_{v-1, n-v-1}, \] so we can obtain a decision rule for testing \(H_0:\left\{\tau_1=\tau_2=\cdots=\tau_v\right\}\) against \(H_A:\left\{\tau_i\right.\) not all equal \(\}\) as reject \(H_0\) if \(\text{msT}/ \text{msE}>F_{v-1, n-v-1, \alpha}\), where \(n-v-1\) is the degrees of freedom for error.
Similarly, the decision rule for testing \(H_0:\{\beta=0\}\) against \(H_A:\{\beta \neq 0\}\), at significance level \(\alpha\), is reject \(H_0\) if \(ms\beta /msE>F_{1, n-v-1, \alpha}\).
\[ \begin{array}{ccccc} \hline \text { Source of variation } & \text { Degrees of freedom } & \text { Sum of squares } & \text { Mean squares } & \text { Ratio } \\ \hline T & v-1 & ssT & \frac{ssT}{v-1} & \frac{msT}{m s E} \\ \beta & 1 & ss\beta & ss\beta & \\ \text { Error } & n-v-1 & s s E & \\ \text { Total } & n-1 & s s_{y y} & m s E \\ \hline \end{array} \]
The inference for a single contrast and multiple contrasts is done similarly to ANOVA. Read Section 9.5 for more details.
The zink experiment was used by C. R. Hicks (1965), Industrial Quality Control, to illustrate the possible bias caused by ignoring an important covariate. The experimental units consisted of 12 steel brackets. Four steel brackets were sent to each of three vendors to be zinc plated. The response variable was the thickness of the zinc plating, in hundred-thousandths of an inch. The thickness of each bracket before plating was measured as a covariate. The data are reproduced in Table 9.8.
data zinc;
input vendor x y;
lines;
1 110 40
1 75 38
1 93 30
1 97 47
2 60 25
2 75 32
2 38 13
2 140 35
3 62 27
3 90 24
3 45 20
3 59 13
;
run;
proc sgplot;
scatter Y=y X=x/group=vendor;
title "Scatter Plot";
run;(B). Fit the ANCOVA model and test the equality of vendor effects.
proc glm data=zinc;
class vendor;
model y=vendor x;
run;\[ \begin{array}{|l|r|r|r|r|r|} \hline \text { Source } & \text { DF } & \text { Type III SS } & \text { Mean Square } & \text { F Value } & \operatorname{Pr}>\text { F } \\ \hline \text { vendor } & 2 & 292.4217278 & 146.2108639 & 3.57 & 0.0778 \\ \hline \mathbf{x} & 1 & 216.1015128 & 216.1015128 & 5.28 & 0.0506 \\ \hline \end{array} \] For \(\alpha=0.05\), we fail to reject the null hypothesis. Therefore there is no significant evidence to indicate a difference among the vendors.
proc glm data=zinc;
class vendor;
model y=vendor;
run;\[ \begin{array}{|r|r|r|r|r|r|} \hline \text { Source } & \text { DF } & \text { Type III SS } & \text { Mean Square } & \text { F Value } & \operatorname{Pr}>\text { F } \\ \hline \text { vendor } & 2 & 665.1666667 & 332.5833333 & 5.51 & 0.0274 \\ \hline \end{array} \]
Now it says there is a significant difference among the vendors! But actually, some of the differences are due to the thickness of the brackets, i.e., the covariate.
proc glm data=zinc;
class vendor;
model y=vendor x vendor*x;
run;\[ \begin{array}{|l|r|r|r|r|r|} \hline \text { Source } & \text { DF } & \text { Type III SS } & \text { Mean Square } & \text { F Value } & \operatorname{Pr}>\text { F } \\ \hline \text { vendor } & 2 & 18.50178803 & 9.25089401 & 0.17 & 0.8450 \\ \hline \mathbf{x} & 1 & 62.68907040 & 62.68907040 & 1.17 & 0.3202 \\ \hline \text { x}^* \text { vendor } & 2 & 6.95209274 & 3.47604637 & 0.07 & 0.9376 \\ \hline \end{array} \] The estimates are provided in the following table. \[ \begin{array}{|l|r|r|r|r|r|} \hline & & & \begin{array}{r} \text { Standard } \\ \text { Error } \end{array} & \text { t Value } & \operatorname{Pr}>|t| \\ \hline \text { Parameter } & \text { Estimate } & & 0.91 & 0.3965 \\ \hline \text { Intercept } & 13.49530957 & \text { B } & 14.78397390 & 0.6464 \\ \hline \text { vendor 1 } & 15.12056598 & \text { B } & 31.31879622 & 0.48 & 0.6464 \\ \hline \text { vendor 2 } & -1.84416489 & \text { B } & 16.98609165 & -0.11 & 0.9171 \\ \hline \text { vendor 3 } & 0.00000000 & \text { B } & & & \\ \hline \text { x } & 0.11726079 & \text { B } & 0.22383267 & 0.52 & 0.6191 \\ \hline \text { x }^* \text { vendor 1 } & -0.00916346 & \text { B } & 0.36785165 & -0.02 & 0.9809 \\ \hline \text { x }^{\star} \text { vendor 2 } & 0.06930605 & \text { B } & 0.24361107 & 0.28 & 0.7856 \\ \hline \text { x }^* \text { vendor 3 } & 0.00000000 & \text { B } & & & \\ \hline \end{array} \] It is interpreted as follows:
The slope estimate for vendor 3 is 0.1172. The slope estimate for vendor 2 is 0.1172+0.0693=0.1865. The slope estimate for vendor 1 is 0.1172-0.0092=0.108.