REPEATED MEASURES ANOVA MATLAB

169 views (last 30 days)
Giulia Soresini
Giulia Soresini on 3 Nov 2019
Commented: Adam Danz on 29 Apr 2021
Hi everyone!
I am having troubles understanding how to perform the correct anova test on my data. I have 16 subjects who are being treated with 2 different treatments. For each subject and each treatment measures are taken at 6 different time points (attached is the excel file). Is repeated measures anova the best way to go? If so, what is the most apporpriate command to use?
Here is the code I have so far. The results are kind of strange so I am not sure if it is wrong or I am just interpreting it the wrong way.
Thanks,
Giulia
Q2 = xlsread('Q2.xls');
Treatment = categorical(Q2(:,2));
ID = Q2(:,1);
t = table(Treatment, Q2(:,3),Q2(:,4),Q2(:,5),Q2(:,6),Q2(:,7), Q2(:,8),'VariableNames', {'Treatment', 't1', 't2', 't3', 't4', 't5', 't6'});
Time = [1 2 3 4 5 6]';
rm = fitrm( t, 't1-t6~Treatment', 'WithinDesign', Time);
ranovatbl = ranova(rm);
disp(ranovatbl);
multcompare(rm, 'Time', 'By', 'Treatment')
  1 Comment
Scott MacKenzie
Scott MacKenzie on 27 Apr 2021
Below is an Anova I just did on your data set (not using MATLAB). What you have is a 2 x 6 within-subjects design. You have two independent variables: Treatment (2 levels) and Time Point (6 levels). What was the dependent variable? In your question, you only say "measures are taken". Measures of what? Here, I'm assuming you measured something like reaction time or task completion time with milliseconds as the units.
Here is a summary of the key results: The effect of treatment on task completion time was not statistically significant, F(1,10) = 1.131, p > .05. The effect of time point on task completion time was statistically significant, F(5,50) = 5.127, p < .001. The Treatment x Time Point interaction effect was also statistically significant, F(5,50) = 6.478, p = .0001. See below for an Anova table for your data and experiment design.
I haven't figure out how to generate a table like this in MATLAB, so I'm using a different tool.
BTW, the attached file has data for 11 participants, not 16 as stated in your question.

Sign in to comment.

Accepted Answer

Adam Danz
Adam Danz on 3 Nov 2019
Edited: Adam Danz on 16 Apr 2021
Is repeated measures anova the best way to go?
Your design is certainly applicable to a repeated measures design but the answer to that question depends on what you're testing.
It looks like you want to perform a repeated measures ANOVA with a covariate. In addition to patient number, the treatment type is the covariate, and the 6 time points is the within-subject repeated measure.
There are 4 null hypotheses to look at:
  • Patient number has no effect on the population mean
  • Treatment type has no effect on the population mean
  • There is no interaction between patient number and treatment type (usually the most important question)
  • There is no 3-way interaction between patient number, treatment type, and treatment time.
Setting up a RM Anova with a convariate in Matlab
Your data are attached in a file named Q2unlocked.xlsx. The following lines set up a RM Anova with a covariate (written in r2019b).
% Read in the table
Q2 = readtable('Q2unlocked.xlsx');
% Replace the 6 measurment headers
Q2.Properties.VariableNames(3:8) = {'t1', 't2', 't3', 't4', 't5', 't6'};
% The within-subjects design
withinDsgn = table((1:6)','VariableNames',{'Time'});
% Run the RM Anova (note the interaction between PATIENT and TREAT!)
rm = fitrm(Q2, 't1-t6~PATIENT*TREAT', 'WithinDesign', withinDsgn);
Also see a similar example provided by Matlab.
Testing for assumptions
Before we look at the results, you must make sure your data are appropriate for a RM Anova by confirming the following assumptions.
  1. Independent observations: You collected the data randomly and without bias
  2. Normality: The measurments have an approximately normal distribution
  3. Sphericity: we'll use the Mauchly test.
% Are the data approximately normal? This should ideally be done for each factor
histogram(reshape(Q2{:,3:end},[],1))
A small rightward tail but not too bad.
% Does the data pass the sphericity test?
rm.mauchly
ans = 1×4 table
W ChiStat DF pValue _______ _______ __ _______ 0.48842 11.537 14 0.64344
Yes: a pValue < 0.05 would fail
RM Anova results
ranova(rm) or rm.ranova produce the results table showing the p-value and then 3 adjusted p-values depending on whether the the response variables have different variances (symmetry assumption).
% Produce ranova table
rm.ranova
ans = 5×8 table
SumSq DF MeanSq F pValue pValueGG pValueHF pValueLB __________ __ __________ ______ ________ ________ ________ ________ (Intercept):Time 3.9654e+06 5 7.9308e+05 2.2798 0.053221 0.071668 0.053221 0.14843 PATIENT:Time 4.0894e+06 5 8.1788e+05 2.3511 0.047012 0.064753 0.047012 0.14259 TREAT:Time 3.2657e+06 5 6.5314e+05 1.8775 0.10606 0.12652 0.10606 0.18747 PATIENT:TREAT:Time 3.9524e+06 5 7.9048e+05 2.2723 0.053916 0.072433 0.053916 0.14905 Error(Time) 3.1309e+07 90 3.4788e+05
Row 1 representing all differences across the within-subjects factors (treatment times). It has a p value of 0.053 which is just beyond the socially accepted threshold of 0.05. A boxplot will show us if this value seems reasonable.
figure()
boxplot(Q2{:,3:end})
xlabel('Treatment times')
ylabel('measurement value')
title('Data pooled between all Patients and treatment factors')
It doesn't look like there's much of an effect of treatment time when factors are combined
Row 2 shows the interaction between Patient and treatment times and is p=0.047
bpdata = [];
for i = 1:max(Q2.PATIENT) %assuming patient numbers are 1:max
bpdata = [bpdata, Q2{Q2.PATIENT==i,3:8},nan(size(unique(Q2.TREAT)))];
end
figure()
boxplot(bpdata)
arrayfun(@xline,7:7:size(bpdata,2))
xlabel('6 treatment times across 11 patients')
ylabel('measurement value')
title('Data pooled between treatment factors')
set(gca,'XTick', [])
Clearly the difference between the 6 treatment times differs across subjects (Treatment type combined)
Row 3 shows the interaction between treatment type and treatment times and is p=0.106.
bpdata = [];
for i = 1:max(Q2.TREAT) %assuming TREAT numbers are 1:max
bpdata = [bpdata, Q2{Q2.TREAT==i,3:8},nan(size(unique(Q2.PATIENT)))];
end
figure()
boxplot(bpdata)
xline(7)
xlabel('6 treatment times across 2 TREAT factors')
ylabel('measurement value')
title('Data pooled between patients')
set(gca,'XTick', [])
Other than the 2nd treatment factor there is no difference between the 2 Treatment types.
Row 4 shows the 3-way interaction between patient, treatment type, and treatment times.
  14 Comments
Adam Danz
Adam Danz on 29 Apr 2021
Thanks for the update, Scott. I saw Jeff's solution earlier today and also noticed the difference in withinDesign. I look forward to examining it closer and perhaps updating my answer when I have a chance to do so, shortly.

Sign in to comment.

More Answers (0)

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!