function [RMAOV1MS] = RMAOV1MS(X,alpha)
%RMAOV1MS Repeated measures single-factor analysis of variance multiple samples
% test.
% This m-file deals with repeated measurements at t time points obtained from
% s groups of subjects or multiple samples. With or without same size. The
% simplest model is:
%
% x_hij = m + g_h +t_j + (gt)_hj + p_i(h) + e_hij
%
% where: x_hij = response at time j from the ith subject in group h; m =
% overall mean; g_h = fixed effect of group h; t_j = fixed effect of
% time j; (gt)_hj = fixed effect for the interaction of the hth group
% with the jth time; p_i(h) = random effects for the ith subject in
% the hth group, and e_hij = independent random error terms.
%
% For we are interested test for differences among groups, test differences
% among time points, and test the significance of the group x time interaction.
% The first test requires the assumption that the within-group covariance
% matrices are equal. The two last-ones tests require this same assumption
% and that the sphericity condition is satisfied. In general, these assumtions
% are required for all tests of within-subjects effects (repeated measures).
%
% Syntax: RMAOV1MS(X,alpha)
%
% Inputs:
% X - data matrix (Size of matrix must be n-by-4;dependent variable=
% column 1;group=column 2;independent variable=column 3;subject=
% column 4).
% alpha - significance level (default = 0.05).
%
% Outputs:
% - Complete Analysis of Variance Table.
%
% Example: From the example 5.4.2 of Charles S. Davis (2002). from an experiment
% comparing dissolution times of pills. We use the logarithms of the
% times, in seconds, for each pill to reach fractions remaining of
% 0.9,0.7,0.5,0.3,0.25, and 0.10. The four groups correspond to
% different storage conditions. Under the assumption thet the
% logarithms of the times are normally distributed, repeated measures
% ANOVA can be used to test for effects due to group, fraction
% remaining, and interaction effect between group and fraction
% remaining. Tha data are shown on the below table.
%
% -----------------------------------------------------------
% Fraction remainig
% ------------------------------------------
% Group Tablet 0.90 0.70 0.50 0.30 0.25 0.10
% -----------------------------------------------------------
% 1 1 2.56 2.77 2.94 3.14 3.18 3.33
% 2 2.64 2.89 3.09 3.26 3.33 3.47
% 3 2.94 3.18 3.33 3.50 3.50 3.66
% 4 2.56 2.83 3.04 3.22 3.26 3.37
% 5 2.64 2.77 2.94 3.14 3.22 3.30
% 6 2.56 2.77 2.94 3.14 3.18 3.26
% 2 1 2.56 2.83 3.07 3.26 3.33 3.51
% 2 2.44 2.74 3.02 3.20 3.28 3.44
% 3 2.34 2.67 2.91 3.16 3.22 3.39
% 4 2.41 2.71 2.94 3.16 3.21 3.36
% 3 1 2.46 2.83 3.09 3.32 3.37 3.54
% 2 2.60 2.93 3.21 3.40 3.46 3.62
% 3 2.48 2.84 3.12 3.35 3.41 3.58
% 4 2.49 2.82 3.05 3.29 3.37 3.52
% 4 1 2.40 2.67 2.94 3.20 3.26 3.47
% 2 2.64 2.94 3.18 3.40 3.45 3.66
% 3 2.40 2.64 2.86 3.09 3.16 3.38
% -----------------------------------------------------------
%
% Data matrix must be:
% X=[2.56 1 1 1;2.77 1 1 2;2.94 1 1 3;3.14 1 1 4;3.18 1 1 5;3.33 1 1 6;2.64 1 2 1;2.89 1 2 2;3.09 1 2 3;3.26 1 2 4;3.33 1 2 5;3.47 1 2 6;
% 2.94 1 3 1;3.18 1 3 2;3.33 1 3 3;3.50 1 3 4;3.50 1 3 5;3.66 1 3 6;2.56 1 4 1;2.83 1 4 2;3.04 1 4 3;3.22 1 4 4;3.26 1 4 5;3.37 1 4 6;
% 2.64 1 5 1;2.77 1 5 2;2.94 1 5 3;3.14 1 5 4;3.22 1 5 5;3.30 1 5 6;2.56 1 6 1;2.77 1 6 2;2.94 1 6 3;3.14 1 6 4;3.18 1 6 5;3.26 1 6 6;
% 2.56 2 1 1;2.83 2 1 2;3.07 2 1 3;3.26 2 1 4;3.33 2 1 5;3.51 2 1 6;2.44 2 2 1;2.74 2 2 2;3.02 2 2 3;3.20 2 2 4;3.28 2 2 5;3.44 2 2 6;
% 2.34 2 3 1;2.67 2 3 2;2.91 2 3 3;3.16 2 3 4;3.22 2 3 5;3.39 2 3 6;2.41 2 4 1;2.71 2 4 2;2.94 2 4 3;3.16 2 4 4;3.21 2 4 5;3.36 2 4 6;
% 2.46 3 1 1;2.83 3 1 2;3.09 3 1 3;3.32 3 1 4;3.37 3 1 5;3.54 3 1 6;2.60 3 2 1;2.93 3 2 2;3.21 3 2 3;3.40 3 2 4;3.46 3 2 5;3.62 3 2 6;
% 2.48 3 3 1;2.84 3 3 2;3.12 3 3 3;3.35 3 3 4;3.41 3 3 5;3.58 3 3 6;2.49 3 4 1;2.82 3 4 2;3.05 3 4 3;3.29 3 4 4;3.37 3 4 5;3.52 3 4 6;
% 2.40 4 1 1;2.67 4 1 2;2.94 4 1 3;3.20 4 1 4;3.26 4 1 5;3.47 4 1 6;2.64 4 2 1;2.94 4 2 2;3.18 4 2 3;3.40 4 2 4;3.45 4 2 5;3.66 4 2 6;
% 2.40 4 3 1;2.64 4 3 2;2.86 4 3 3;3.09 4 3 4;3.16 4 3 5;3.38 4 3 6];
%
% Calling on Matlab the function:
% RMAOV1MS(X)
%
% Answer is:
%
% The number of groups are: 4
% The maximum number of subjects in a group are: 6
% The minimum number of subjects in a group are: 3
% The number of IV levels are: 6
%
% ANOVA table for repeated measures single-factor multiple samples.
% ------------------------------------------------------------------------------
% SOV SS df MS F P
% ------------------------------------------------------------------------------
% Groups 0.2038 3 0.0679 0.87 0.4798
% Subjts(Gps.) 1.0108 13 0.0778
% IV 10.0745 5 2.0149 3198.78 0.0000
% Group x IV 0.2045 15 0.0136 21.64 0.0000
% Residual 0.0409 65 0.0006
% ------------------------------------------------------------------------------
% Total 11.5344 101
% ------------------------------------------------------------------------------
% If the P results are smaller than 0.05
% the corresponding Ho's tested result statistically significant. Otherwise,
% are not significative.
%
% Created by A. Trujillo-Ortiz, R. Hernandez-Walls, K. Barba-Rojo and
% A. Castro-Perez
% Facultad de Ciencias Marinas
% Universidad Autonoma de Baja California
% Apdo. Postal 453
% Ensenada, Baja California
% Mexico.
% atrujo@uabc.mx
%
% Copyright.December 13, 2006.
%
% To cite this file, this would be an appropriate format:
% Trujillo-Ortiz, A., R. Hernandez-Walls, K. Barba-Rojo and A. Castro-Perez. (2006).
% RMAOV1MS:Repeated measures single-factor analysis of variance multiple samples test.
% A MATLAB file. [WWW document]. URL http://www.mathworks.com/matlabcentral/
% fileexchange/loadFile.do?objectId=13380
%
% Reference:
% Davis, S. Ch. (2002), Statistical Methods for the Analysis of Repeated
% Measurements. Springer-Verlag:New York. Chapter 5.
%
if nargin < 2,
alpha = 0.05; %(default)
end;
if (alpha <= 0 | alpha >= 1)
fprintf('Warning: significance level must be between 0 and 1\n');
return;
end;
if nargin < 1,
error('Requires at least one input argument.');
return;
end;
a = max(X(:,2)); %Number of groups.
b = max(X(:,3)); %Maximum number of subjects in a group.
t = max(X(:,4)); %Number of levels of independent variable.
indice = X(:,2);
for i = 1:a
Xe = find(indice==i);
eval(['G' num2str(i) ' = X(Xe,1);']);
end
indice = X(:,4);
for k = 1:t
Xe = find(indice==k);
eval(['IV' num2str(k) '=X(Xe,1);']);
end
C = (sum(X(:,1)))^2/length(X(:,1)); %Correction term.
SSTO = sum(X(:,1).^2)-C; %Total sum of squares.
dfTO = length(X(:,1))-1; %Total degrees of freedom.
for i = 1:a
for j = 1:b
Xe = find((X(:,2)==i) & (X(:,3)==j));
eval(['T' num2str(i) num2str(j) ' = X(Xe,1);']);
end
end
for i = 1:a
for k = 1:t
Xe = find((X(:,2)==i) & (X(:,4)==k));
eval(['GI' num2str(i) num2str(k) '=X(Xe,1);']);
end
end
%Procedure related to the subjects in groups.
NSG = [];
for i = 1:a
for j = 1:b
eval(['n' num2str(i) num2str(j) ' = length(T' num2str(i) num2str(j) ');'])
eval(['x = n' num2str(i) num2str(j) ';'])
NSG = [NSG,x];
end
end
%Due it is an unequal sample sizes design, some subgroups arithmetic operations
%could been not-a-number (like 0.0/0.0). So it is necessary to make the estimations
%ignoring NaNs.
warning off;
SG = [];
for i = 1:a
for j = 1:b
eval(['x = ((sum(T' num2str(i) num2str(j) ').^2)/length(T' num2str(i) num2str(j) '));'])
SG = [SG,x];
nans = isnan(SG);
r = find(nans);
SG(r) = zeros(size(r));
end
end
warning on;
SSSG = sum(SG)-C;
dfSG = (a*b)-1;
%Procedure related to the groups.
G = [];
for i = 1:a
eval(['x = ((sum(G' num2str(i) ').^2)/length(G' num2str(i) '));'])
G = [G,x];
end
SSG = sum(G)-C;
dfG = a-1;
MSG = SSG/dfG;
%Procedure related to the IV (independent variable).
IV = [];
for k = 1:t
eval(['x =((sum(IV' num2str(k) ').^2)/length(IV' num2str(k) '));']);
IV = [IV,x];
end
SSIV = sum(IV)-C; %sum of squares for the IV
dfIV = t-1; %degrees of freedom for the IV
MSIV = SSIV/dfIV; %mean square for the IV
%Procedure related to the interaction group-independent variable.
GI = [];
for i = 1:a
for k = 1:t
eval(['x =((sum(GI' num2str(i) num2str(k) ').^2)/length(GI' num2str(i) num2str(k) '));']);
GI = [GI,x];
end
end
SSGI = sum(GI)-sum(G)-sum(IV)+C; %sum of squares of the Group x IV
dfGI = dfG*dfIV; %degrees of freedom of the Group x IV
MSGI = SSGI/dfGI; %mean square for the Group x IV
%Procedure for estimation of sample size for unequal sample sizes.
no1 = []; no2 = [];
for i = 1:a
for j = 1:b
eval(['x = length(T' num2str(i) num2str(j) ');']);
no1 = [no1,x];
eval(['x = length(T' num2str(i) num2str(j) ').^2;']);
no2 = [no2,x];
end
end
%Procedure related to the within subgroups (error).
SSE = SSTO-SSSG;
nneg = length(find((no1-1)<0));
dfE = sum(no1-1)+nneg;
MSE = SSE/dfE;
%Procedure related to the among subgroups within groups.
SSSGG = SSTO-SSG-SSE;
dfSGG = dfTO-dfG-dfE;
MSSGG = SSSGG/dfSGG;
%Procedure related to the residual.
SSR = SSTO-SSG-SSSGG-SSIV-SSGI;
dfR = dfTO-dfG-dfSGG-dfIV-dfGI;
MSR = SSR/dfR;
fprintf('The number of groups are:%2i\n', a);
fprintf('The maximum number of subjects in a group are:%2i\n', b);
c = min(full(sum(spones(reshape(no2,b,a))))); %Minimum number of subjects in a group.
fprintf('The minimum number of subjects in a group are:%2i\n', c);
fprintf('The number of IV levels are:%2i\n\n', t);
FG = MSG/MSSGG;
FIV = MSIV/MSR;
FGI = MSGI/MSR;
PrG = 1 - fcdf(FG,dfG,dfSGG); %Probability associated to the groups F statistic.
PrIV = 1 - fcdf(FIV,dfIV,dfR); %Probability associated to the IV F statistic.
PrGI = 1 - fcdf(FGI,dfGI,dfR); %Probability associated to the interaction group-IV F statistic.
disp('ANOVA table for repeated measures single-factor multiple samples.')
fprintf('------------------------------------------------------------------------------\n');
disp('SOV SS df MS F P ')
fprintf('------------------------------------------------------------------------------\n');
fprintf('Groups %18.4f%10i%14.4f%13.2f%9.4f\n',SSG,dfG,MSG,FG,PrG);
fprintf('Subjts(Gps.)%19.4f%10i%14.4f\n',SSSGG,dfSGG,MSSGG);
fprintf('IV %25.4f%10i%14.4f%13.2f%9.4f\n',SSIV,dfIV,MSIV,FIV,PrIV);
fprintf('Group x IV %20.4f%10i%14.4f%13.2f%9.4f\n',SSGI,dfGI,MSGI,FGI,PrGI);
fprintf('Residual %22.4f%10i%14.4f\n',SSR,dfR,MSR);
fprintf('------------------------------------------------------------------------------\n');
fprintf('Total %25.4f%10i\n',SSTO,dfTO);
fprintf('------------------------------------------------------------------------------\n');
fprintf('If the P results are smaller than% 3.2f\n', alpha );
disp('the corresponding Ho''s tested result statistically significant. Otherwise,');
disp('are not significative.');
disp(' ')
return,