ANOVA WITH REPEATED MEASURES FACTORS
Preliminaries
Model Formulae
If you have not yet read the Model Formulae tutorial, you should do so now.
Some Explanations
The naming of experimental designs has led to a great deal of confusion on the part of students as well as tutorial writers. The first design I will discuss in this tutorial is the single factor within subjects design, also called the single factor repeated measures design. This is a design with one response variable, and each "subject" or "experimental unit" is measured multiple times on that response variable. For example, pre-post designs, in which a group of subjects is measured both before and after some treatment, fit into this category. This design is also referred to in some sources as a "treatment × subjects" design (read "treatment by subjects"). You can also think of the design as having treatments nested within subjects. Or you could also think of it as a block design in which subjects correspond to blocks. You see my point.
In my field (psychology), we are used to referring to our experimental units as "subjects," because they are almost always either human beings or animals. I will retain that terminology in this tutorial. This is why, in the first study discussed below, grocery items are referred to as subjects--because the groceries correspond to what would ordinarily be human subjects in a typical psychology experiment. For a change, it is people in other fields who will have to adjust to strange terminology!
Single Factor Designs With Repeated Measures
Several years ago a friend and I decided to compare prices in local grocery
stores. We made a list of ten representative grocery items and then went to four
local stores and recorded the price on each item in each store. The most
convenient way to table the data--and the first way that comes to mind--is as
follows (numbers are prices in dollars)...
subject storeA storeB storeC storeD
lettuce 1.17 1.78 1.29 1.29
potatoes 1.77 1.98 1.99 1.99
milk 1.49 1.69 1.79 1.59
eggs 0.65 0.99 0.69 1.09
bread 1.58 1.70 1.89 1.89
cereal 3.13 3.15 2.99 3.09
ground.beef 2.09 1.88 2.09 2.49
tomato.soup 0.62 0.65 0.65 0.69
laundry.detergent 5.89 5.99 5.99 6.99
aspirin 4.46 4.84 4.99 5.15
Excuse for the moment the absence of nicities in this table, like stubs and
spanners, because shortly you will be happy I did the table this way!
First, I remember when a can of tomato soup was 12¢, but that's another story. We should notice laundry detergent costs a lot more than soup. This isn't interesting to us. In other words, there is a great deal of variability--perhaps 99% or more of the total variability--attributable entirely to individual differences between subjects, in which we are not the least bit interested. If we treated this as a straightforward oneway randomized ANOVA, all of this variability would go into the error term, and any variability due to store would be totally swamped. The advantage of the repeated measures analysis is that it allows us to parcel out variability due to subjects.
Second, we are now faced with getting these data into R. Try this...
> url = "http://ww2.coastal.edu/kingw/statistics/R-tutorials/text/groceries.txt"
> groceries = read.table(url, header=T)
> groceries
Whether or not that will work will depend upon whose servers are up and how much
access your server is granted to coastal.edu. Ordinarily, you should be able to
read data tables from a website, but who knows how these things sometimes work?
The facilities to do so are provided in R.
You can also try it this way. Click on this
link, and when the data table comes up in your browser, copy and paste it
to a text editor, save the file as "groceries.txt", drop it into your working
directory, and then...
> groceries = read.table("groceries.txt", header=T)
> groceries
Or, you can go back to that link, copy the entire table, and then paste it at
the "0:" prompt...
> groceries = read.table(stdin(), header=T)
0: subjectstoreAstoreBstoreCstoreD
1: lettuce1.171.781.291.29
2: potatoes1.771.981.991.99
3: milk1.491.691.791.59
4: eggs0.650.990.691.09
5: bread1.581.701.891.89
6: cereal3.133.152.993.09
7: ground.beef2.091.882.092.49
8: tomato.soup0.620.650.650.69
9: laundry.detergent5.895.995.996.99
10: aspirin4.464.844.995.15
11:
> groceries
subjectstoreAstoreBstoreCstoreD
1 lettuce1.171.781.291.29
2 potatoes1.771.981.991.99
3 milk1.491.691.791.59
4 eggs0.650.990.691.09
5 bread1.581.701.891.89
6 cereal3.133.152.993.09
7 ground.beef2.091.882.092.49
8 tomato.soup0.620.650.650.69
9 laundry.detergent5.895.995.996.99
10 aspirin4.464.844.995.15
As you can see, that didn't work for me because R didn't pick up the tabs as
being white space. What can I tell ya? Sometimes that works, and sometimes it
doesn't, and I don't know why! More R idiosyncrasies! (Hey! It's free!) If the
same thing happened to you, try this. Copy and paste this script into R
(WARNING: this will erase all items in your working directory)...
### start copying with this line
rm(list=ls())
groceries = data.frame(
c("lettuce","potatoes","milk","eggs","bread","cereal",
"ground.beef","tomato.soup","laundry.detergent","aspirin"),
c(1.17,1.77,1.49,0.65,1.58,3.13,2.09,0.62,5.89,4.46),
c(1.78,1.98,1.69,0.99,1.70,3.15,1.88,0.65,5.99,4.84),
c(1.29,1.99,1.79,0.69,1.89,2.99,2.09,0.65,5.99,4.99),
c(1.29,1.99,1.59,1.09,1.89,3.09,2.49,0.69,6.99,5.15))
colnames(groceries) = c("subject","storeA","storeB","storeC","storeD")
groceries
### end copying with this line
And if that didn't do it, well, start typing!
This is not a proper data frame. Why not? Because our response variable--our
ONE response variable--is the price of the items, and this is listed in four
columns. Very bad! We need each variable in ONE column of the data frame. This
sort of problem occurs so often, you'd think some nice R coder would have
written a utility to rearrange such tables. (Remember the pleasure we had
rearranging the "anorexia" data frame in an earlier tutorial?) Well, almost...
> stack(groceries) # Funny! :^)
I have not shown the output, but this is ALMOST what we want. It would have
worked beautifully if this were a between subjects design. Here it is but a
start...
> groceries2 = stack(groceries)
> subject = rep(groceries$subject,4) # create the "subject" variable
> groceries2[3] = subject # add it to the new data frame
> rm(subject) # clean up your workspace
> colnames(groceries2) = c("price", "store", "subject") # rename the columns
> groceries2 # take a look
price store subject
1 1.17 storeA lettuce
2 1.77 storeA potatoes
3 1.49 storeA milk
4 0.65 storeA eggs
5 1.58 storeA bread
6 3.13 storeA cereal
7 2.09 storeA ground.beef
8 0.62 storeA tomato.soup
9 5.89 storeA laundry.detergent
10 4.46 storeA aspirin
11 1.78 storeB lettuce
12 1.98 storeB potatoes
13 1.69 storeB milk
14 0.99 storeB eggs
15 1.70 storeB bread
16 3.15 storeB cereal
17 1.88 storeB ground.beef
18 0.65 storeB tomato.soup
19 5.99 storeB laundry.detergent
20 4.84 storeB aspirin
21 1.29 storeC lettuce
22 1.99 storeC potatoes
23 1.79 storeC milk
24 0.69 storeC eggs
25 1.89 storeC bread
26 2.99 storeC cereal
27 2.09 storeC ground.beef
28 0.65 storeC tomato.soup
29 5.99 storeC laundry.detergent
30 4.99 storeC aspirin
31 1.29 storeD lettuce
32 1.99 storeD potatoes
33 1.59 storeD milk
34 1.09 storeD eggs
35 1.89 storeD bread
36 3.09 storeD cereal
37 2.49 storeD ground.beef
38 0.69 storeD tomato.soup
39 6.99 storeD laundry.detergent
40 5.15 storeD aspirin
NOW we have a proper data frame for this analysis. I belabor this point because
it is a critical one!
At which store should we shop?
> with(groceries2, tapply(price, store, sum))
storeA storeB storeC storeD
22.85 24.65 24.36 26.26
For the items in the sample, it looks like storeA had the lowest prices. But as
this is only a random sample of items we might have wanted to buy, we have to
ask, will this difference hold up in general?
Once the data frame is set up correctly (all values of the response variable
in ONE column, and another column holding the subject variable--i.e., which
identifies which values of the response go with which subjects), the trick to
doing a repeated measures ANOVA is in getting the model formula right. The
Error term must reflect that we have "treatments nested within subjects." That
is to say, in principle at least, we can see the "store" effect within each and
every "subject" (grocery item). The ANOVA is properly done this way...
> aov.out = aov(price ~ store + Error(subject/store), data=groceries2)
> summary(aov.out)
Error: subject
Df Sum Sq Mean Sq F value Pr(>F)
Residuals 9 115.193 12.799
Error: subject:store
Df Sum Sq Mean Sq F value Pr(>F)
store 3 0.58586 0.19529 4.3442 0.01273 *
Residuals 27 1.21374 0.04495
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
"Store" and "subject" are our sources of variability. The treatment we are
interested in is "store" (that's what we want to see the effect of), and this
treatment effect is visible within each subject (i.e., nested within each
subject). So the proper Error term is "subject/store", which is read as "store
within subject." Notice once all that "subject" variability is parceled out, we
do have a significant "store" main effect.
The Tukey HSD test will not run on this model object, so I propose pairwise
t-tests with adjusted p-values to see which stores are significantly
different (note: the output will be a table of p-values for each possible
pairwise comparison)...
> with(groceries2, pairwise.t.test(price, store,
+ p.adjust.method="holm", paired=T))
Pairwise comparisons using paired t tests
data: price and store
storeA storeB storeC
storeB 0.17 - -
storeC 0.17 0.69 -
storeD 0.07 0.49 0.33
P value adjustment method: holm
The syntax for pairwise.t.test( ) is a bit out of line with other R
functions, so this is how it works. The first argument is the response variable.
The second argument is the grouping or explanatory variable. There is a comma
between these two arguments--the function does not accept a formula. There are
many p.adjust.methods available. See ?p.adjust for details. Since these should
be paired (i.e., dependent) t-tests, I set the "paired=" option to TRUE. We see
no significant differences in our post hoc tests, but storeA vs. storeD was at
least close, p = .07.
The Friedman Rank Sum Test
There is a nonparametric version of the oneway ANOVA with repeated measures.
It is executed this way...
> friedman.test(price ~ store | subject, data=groceries2)
Friedman rank sum test
data: price and store and subject
Friedman chi-squared = 13.3723, df = 3, p-value = 0.003897
The formula reads "price as a function of store given subject". Here, "subject"
is being treated as a blocking variable.
Two Way Mixed Factorial Designs
Clean out your workspace, and then enter the data for this analysis by
copying and pasting this script into your R Console...
### begin copying with this line
# create the data vectors for the four conditions
pre.chicken=c(18,13,18,15,22,32,31,24,15)
post.chicken=c(15,13,17,15,24,31,31,25,17)
pre.pasta=c(17,30,18,13,23,27,27,24,23)
post.pasta=c(19,31,18,13,24,27,26,28,26)
# create the response variable vector
SSS=c(pre.chicken,post.chicken,pre.pasta,post.pasta)
# create the test condition labels
test=rep(c("pre","post","pre","post"),c(9,9,9,9))
# create the diet condition labels
diet=rep(c("chicken","pasta"),c(18,18))
# create the subject labels
subj.chicken=rep(c("A","B","C","D","E","F","G","H","I"),2)
subj.pasta=rep(c("J","K","L","M","N","O","P","Q","R"),2)
# create the subject vector
subject=c(subj.chicken,subj.pasta)
# create the dataframe
hill=data.frame(subject,diet,test,SSS)
# clean up
rm(pre.chicken,post.chicken,pre.pasta,post.pasta)
rm(SSS,test,diet,subj.chicken,subj.pasta,subject)
### end copying with this line
You should now have a data frame called "hill" in your workspace...
> ls()
[1] "groceries" "groceries2" "hill"
> str(hill)
'data.frame': 36 obs. of 4 variables:
$ subject: Factor w/ 18 levels "A","B","C","D",..: 1 2 3 4 5 6 7 8 9 1 ...
$ diet : Factor w/ 2 levels "chicken","pasta": 1 1 1 1 1 1 1 1 1 1 ...
$ test : Factor w/ 2 levels "post","pre": 2 2 2 2 2 2 2 2 2 1 ...
$ SSS : num 18 13 18 15 22 32 31 24 15 15 ...
These data were collected by Ben Hill (by now, Ben Hill, Ph.D.) as part of his
senior research project conducted in our department during the Fall 1999
semester. His interest was in the role of brain serotonin in creating sensation
seeking behavior. Since he couldn't manipulate brain serotonin directly (our
IRB would not permit it!), he chose to do so by way of diet, relying on the
finding that brain serotonin is elevated after a meal rich in carbohydrates.
Therefore, he had 18 subjects come to the lab. All subjects were given the
Sensation Seeking Survey. Half the subjects were then given an evening meal
consisting primarily of chicken, while the other half were given a meal
consisting primarily of pasta. An hour later, the Sensation Seeking Survey was
readministered. In the data frame, "subject" identifies the subject so that
his or her "pre" and "post" eating scores can be kept paired for the analysis,
"diet" identifies the type of meal the subject was given, "test" contains the
pre-post meal information, and SSS contains the survey scores.
In this experiment, we have two explanatory variables, "diet", which is a
between-subjects variable, and "test", which is a within-subjects variable.
Thus, in lingo familiar to social scientists, the design is called a
2 × 2 mixed factorial design with repeated measures on "test".
It's worth making sure you understand how the data frame is structured, so let's
have a look...
> hill
subject diet test SSS
1 A chicken pre 18
2 B chicken pre 13
3 C chicken pre 18
4 D chicken pre 15
5 E chicken pre 22
6 F chicken pre 32
7 G chicken pre 31
8 H chicken pre 24
9 I chicken pre 15
10 A chicken post 15
11 B chicken post 13
12 C chicken post 17
13 D chicken post 15
14 E chicken post 24
15 F chicken post 31
16 G chicken post 31
17 H chicken post 25
18 I chicken post 17
19 J pasta pre 17
20 K pasta pre 30
21 L pasta pre 18
22 M pasta pre 13
23 N pasta pre 23
24 O pasta pre 27
25 P pasta pre 27
26 Q pasta pre 24
27 R pasta pre 23
28 J pasta post 19
29 K pasta post 31
30 L pasta post 18
31 M pasta post 13
32 N pasta post 24
33 O pasta post 27
34 P pasta post 26
35 Q pasta post 28
36 R pasta post 26
Note that all the SSS scores are in ONE column! The other columns completely
identify the experimental conditions associated with that score, including
which subject it came from.
By way of data summary...
> with(hill, tapply(SSS, list(diet,test), mean))
post pre
chicken 20.88889 20.88889
pasta 23.55556 22.44444
Drat! There goes R arranging our factor levels in alphabetical order again,
making our means table look backwards. Let's pretty that up...
> hill$test = factor(hill$test, levels=c("pre","post"))
> with(hill, tapply(SSS, list(diet,test), mean))
pre post
chicken 20.88889 20.88889
pasta 22.44444 23.55556
Better! A nice graph might be in order as well...
> with(hill, boxplot(SSS ~ diet + test)) # output not shown
> with(hill, boxplot(SSS ~ test + diet)) # compare this one with the last one
> title(main="Ben Hill's SSS Data")
> title(ylab="SSS Scores")
Once again, the trick is in getting the model formula correct. In this case,
we have two explanatory variables, and we want to see all possible main effects
and interactions. The "test" variable is "within subjects"...
> aov.out = aov(SSS ~ diet * test + Error(subject/test), data=hill)
> summary(aov.out)
Error: subject
Df Sum Sq Mean Sq F value Pr(>F)
diet 1 40.11 40.11 0.5094 0.4857
Residuals 16 1259.78 78.74
Error: subject:test
Df Sum Sq Mean Sq F value Pr(>F)
test 1 2.7778 2.7778 2.1739 0.1598
diet:test 1 2.7778 2.7778 2.1739 0.1598
Residuals 16 20.4444 1.2778
The effect we're looking for would be shown by the diet × test
interaction. And sadly, it's not significant.
Other Designs
For reference, here are model formulae for a couple other common designs...
Two factor design with repeated measures on both factors:
DV ~ IV1 * IV2 + Error(subject/(IV1*IV2))
Three factor mixed design with repeated measures on IV2 and IV3:
DV ~ IV1 * IV2 * IV3 + Error(subject/(IV2*IV3))
And so on.