Chapter 32 Interaction Plots

32.1 Set up

We’re going to look at lots of examples of interaction plots. These plots are just “connect the dots” plots of tables of means, so we’ll start with a bit about getting those tables.

We’ll begin looking at the SproutingBarley data set from Table 8.1. Two factors (amount of water and age of seeds); response is number of seeds sprouting.

> library(cfcdae)
> data(SproutingBarley)
> head(SproutingBarley)
  weeks water sprouting weeks.z water.z
1     1     4        11       1       4
2     1     8         8       1       8
3     3     4         7       3       4
4     3     8         1       3       8
5     6     4         9       6       4
6     6     8         5       6       8

The function tapply() allows us to get summaries of its first argument by the levels of the factors listed in its second argument. The result is a vector (for one factor), a matrix (for two factors), and so on with dimensions corresponding to the numbers of levels of the factors in the list. The third argument is the function you want applied to the groups after the data have been split, in this case, the mean. These are the means that we will plot in interaction plots.

> means1 <- with(SproutingBarley,tapply(sprouting,list(weeks,water),mean))
> means1
           4         8
1   8.666667  4.666667
3  13.333333  3.666667
6  21.000000  7.666667
9  25.333333  6.666667
12 34.000000 17.000000

We can do it the other way too.

> means2 <- with(SproutingBarley,tapply(sprouting,list(water,weeks),mean))
> means2
         1         3         6         9 12
4 8.666667 13.333333 21.000000 25.333333 34
8 4.666667  3.666667  7.666667  6.666667 17

Of course, one is just the transpose of the other.

> t(means1)
         1         3         6         9 12
4 8.666667 13.333333 21.000000 25.333333 34
8 4.666667  3.666667  7.666667  6.666667 17

You can get interaction plots using the interactplot() function in cfcdae. Base R includes an interaction.plot() function, and the gplots package includes a plotmeans() function. interactplot() is based on interaction.plot() but contains some of the features of plotmeans().

The way it works is that the first argument is the factor used along the horizontal axis, the second argument defines the traces that get drawn, and the third argument gives the response to be averaged. There are a lot more optional arguments, some of which we will eventually use.

32.2 Sprouting barley plots

There appears to be some interaction: as we move to older seeds, the gaps between the levels of water tend to get bigger.

> interactplot(weeks,water,sprouting,data=SproutingBarley)
Base interaction plot of sprouting barley

Figure 32.1: Base interaction plot of sprouting barley

Interaction plots can look very differently depending on which variable is put on the x-axis and which is used to label the lines. Here we reverse the order. I find this plot less understandable than the first one, even though it contains exactly the same information. In this plot, the slopes get steeper as we move to higher levels of weeks; that’s like the gap getting bigger in the previous plot.

> interactplot(water,weeks,sprouting,data=SproutingBarley)
Transposed interaction plot of sprouting barley

Figure 32.2: Transposed interaction plot of sprouting barley

In general, I tend to prefer interaction plots with fewer lines and more horizontal levels.

You can also add confidence intervals for each mean. By default we use a pooled estimate of error computed by pooling the variances within each of the groups (essentially the MSE in the separate means model).

In this example, the confidence intervals between levels of water go from mostly overlapping to mostly not overlapping, so the interaction is probably real and not just random noise.

> interactplot(weeks,water,sprouting,confidence=.95,data=SproutingBarley)
Interaction plot with confidence intervals for Barley data

Figure 32.3: Interaction plot with confidence intervals for Barley data

If you use pooled=FALSE, you get confidence intervals based on the variance in each group. These estimates of variance will have fewer degrees of freedom, and the t-intervals could be much wider because of that.

Here we can see signs of nonconstant variance.

> interactplot(weeks,water,sprouting,confidence=.95,pooled=FALSE,data=SproutingBarley)
Interaction plot with confidence intervals for barley data using individual treatment variances

Figure 32.4: Interaction plot with confidence intervals for barley data using individual treatment variances

If we transform the data, the nonconstant variance problem is reduced, and the interaction is not as pronounced.

> interactplot(weeks,water,sqrt(sprouting),confidence=.95,pooled=FALSE,data=SproutingBarley)
Interaction plot with confidence intervals for square root barley data using individual treatment variances

Figure 32.5: Interaction plot with confidence intervals for square root barley data using individual treatment variances

32.3 Chick body weight

Chick body weights as a function of phosphorus and calcium factors in their food. It’s already a table of means. If there were replications in the original data, we don’t have them.

> data(ChickBodyWeight)
> ChickBodyWeight
  weight    P  Ca
1    584 0.25 0.6
2    489 0.25 0.9
3    453 0.25 1.2
4    616 0.50 0.6
5    606 0.50 0.9
6    621 0.50 1.2

For phosphorus 1, means change with levels of calcium, but they are fairly steady for phosphorus 2. Again, this is all visual with no inference behind it.

> interactplot(Ca,P,weight,data=ChickBodyWeight)
Warning in interactplot(Ca, P, weight, data = ChickBodyWeight): no degrees of freedom for estimating
error
Interaction for chick body weight; Ca horizontal

Figure 32.6: Interaction for chick body weight; Ca horizontal

I think the first version is easier to understand that putting phosphorus on the horizontal.

> interactplot(P,Ca,weight,data=ChickBodyWeight)
Warning in interactplot(P, Ca, weight, data = ChickBodyWeight): no degrees of freedom for estimating
error
Interaction for chick body weight; P horizontal

Figure 32.7: Interaction for chick body weight; P horizontal

32.4 Zinc retention

This is the Hunt and Larson liver zinc data. This will be our first more-than-two-way example. The treatment factors are zinc in the diet, zinc in the last meal, and protein in the last meal. Again, these data are either treatment means or the data from a single replication.

> data(ZincRetention)
> head(ZincRetention)
  m.zinc d.zinc m.protein retention
1      2      1         1        48
2      2      1         2        44
3      2      1         3        59
4      2      1         4        64
5      1      1         1        52
6      1      1         2        73
> interactplot(m.protein,d.zinc,retention,data=ZincRetention)
Zinc retention interaction; protein horizontal, diet zinc trace

Figure 32.8: Zinc retention interaction; protein horizontal, diet zinc trace

Warning: if you ask for confidence intervals for a plot that does not include all of the factors in the experiment, the intervals will usually be too long. This is because ignoring some factors will lead to using an MSE in the confidence intervals that is too large (all the variability due to the ignored factors is inflating the MSE).

> interactplot(m.protein,d.zinc,retention,data=ZincRetention,
+     confidence=.95,main="These intervals are too long")
Zinc retention interaction with incorrect CI; protein horizontal, diet zinc trace

Figure 32.9: Zinc retention interaction with incorrect CI; protein horizontal, diet zinc trace

You can fix this by telling interactplot() what to use for the MSE and how many df it is based on. In a few lectures you’ll be able to figure out where the magically appearing sigma2 and df actually came from. In any case, the evidence for interaction is fairly weak.

> interactplot(m.protein,d.zinc,retention,data=ZincRetention,
+   confidence=.95,sigma2=9,df=7)
Zinc retention interaction with CI; protein horizontal, diet zinc trace

Figure 32.10: Zinc retention interaction with CI; protein horizontal, diet zinc trace

The meal zinc by diet zinc seems almost perfectly additive.

> interactplot(m.zinc,d.zinc,retention,data=ZincRetention)
Zinc retention interaction with CI; meal zinc horizontal, diet zinc trace

Figure 32.11: Zinc retention interaction with CI; meal zinc horizontal, diet zinc trace

Meal protein by meal zinc shows interaction. I have added the intervals to remove doubt about it being just noise.

> interactplot(m.protein,m.zinc,retention,data=ZincRetention,
+   confidence=.95,sigma2=9,df=7)
Zinc retention interaction with CI; protein horizontal, meal zinc trace

Figure 32.12: Zinc retention interaction with CI; protein horizontal, meal zinc trace

You can combine factors to make the trace.

> interactplot(m.protein,m.zinc:d.zinc,retention,data=ZincRetention)
Warning in interactplot(m.protein, m.zinc:d.zinc, retention, data = ZincRetention): no degrees of
freedom for estimating error
Zinc retention interaction; protein horizontal, meal zinc by diet zinc trace

Figure 32.13: Zinc retention interaction; protein horizontal, meal zinc by diet zinc trace

You can combine factors to make the horizontal

> interactplot(m.zinc:d.zinc,m.protein,retention,data=ZincRetention)
Warning in interactplot(m.zinc:d.zinc, m.protein, retention, data = ZincRetention): no degrees of
freedom for estimating error
Zinc retention interaction; meal zinc by diet zinc horizontal, protein trace

Figure 32.14: Zinc retention interaction; meal zinc by diet zinc horizontal, protein trace