Convert anovan syntax to linear mixed-effect model - Nested experimental design

I would like to convert the syntax below using anovan into a linear mixed-effect model. Below a mock and small example of my experimental design. Namely, I have two genotypes (WT and Mut), for each genotype I choose randomly a mouse (Mouse ID), and from each mouse I take slices (SliceID), and finally from each slice I obtain a slope value (Slope). To summarize, I have a fixed factor (Genotype), random variable MouseID nested in Genotype and random variable SliceID nested in MouseID variable. My experimental question is the following: Does slope (my dependent variable) change across genotypes? It is a nested design experiment, and I would like to reflect the desing in my linear mixed-effect model.

Slope=[12 15 16 10 100 102 135 107];
SliceID=[{'s1'},{'s2'},{'s3'},{'s4'},{'s13'},{'s14'},{'s15'},{'s16'}];
MouseID=[{'M1'},{'M1'},{'M2'},{'M2'},{'M3'},{'M3'},{'M4'},{'M4'}];
Genotype=[{'WT'},{'WT'},{'WT'},{'WT'},{'Mut'},{'Mut'},{'Mut'},{'Mut'}];
nestedMatrix=[0 0 0; 1 0 0; 0 1 0];
[p,tbl,stats,terms] = anovan(Slope,{Genotype,MouseID, SliceID},'random', [2,3],'nested',nestedMatrix,'varnames',{'Genotype','MouseID','SliceID'});

DataTable=table(Genotype',MouseID', SliceID',Slope','VariableNames',{'Genotype','MouseID','SliceID','Slope'});

I would like to convert the anovan syntax into a linear mixed-effect model (using the function fitlme). This is because the anovan function is not supported for further multicompare in case I would like also to introduce a repetitive measure (longitudinal study). My difficulty is to code for the fitlme function the appropriate syntax to reflect two random and nested variables (SliceID nested in MouseID and MouseID nested in Genotype).

I tried to use the examples in the Mathwork pages, but none really fit my design.

Any help is much appreciated

Answer 1

[Mattia Bonzanni]

Iwas able to have the same result using both syntaxes.

With the mock example above, even if the F-statistic is correct, the anova on the lme model gave no pValues given df2=0 using the 'DFMethod','satterthwaite' syntax.

I attached a larger mock example (DatasetTry2.mat; G:Genotype; M:MouseID; S:SliceID) and I compared the results from the anovan with the experimental design described above and the lme function as follows:

load DatasetTry2
nestedMatrix=[0 0 0; 1 0 0; 0 1 0]; % nested matrix for the anovan function
[p,tbl,stats,terms] = anovan(DatasetTry{:,4},{DatasetTry{:,1},DatasetTry{:,2}, DatasetTry{:,3}},'random', [2,3],'nested',nestedMatrix,'varnames',{'Genotype','MouseID','SliceID'});

lme=fitlme(DatasetTry,'Slope~G+(1|M)','FitMethod','REML'); % I used the Fit method REML(Restricted maximum likelihood estimation) considering that REML is the default in most cases because it provides unbiased estimates, while ML is only unbiased for large designs
statsMM2 = anova(lme,'DFMethod','satterthwaite'); % I used the Satterthwaite approximation since it can be used if the variances are either equal or different. 

From the anovan syntax, the Genotype fixed factor had F=193.82 and a pVal=7.15e-0.8. I realized that it is not necessary for my experimental question, namely if the slope values are equal between genotypes (WT and APP in the dataset), to decleare the SliceID variable. Indeed, the follwoing code gave the same F and p values for the Genotype factor:

nestedMatrix=[0 0 ; 1 0];
[p,tbl,stats,terms] = anovan(DatasetTry{:,4},{DatasetTry{:,1},DatasetTry{:,2}},'random', [2],'nested',nestedMatrix,'varnames',{'Genotype','MouseID'}); % omitting the random SliceID variable

Thus, I did not use the random variable Slice (S) in the design of the linear mixed-effect model. I was thus left with Genotype (G: WT and APP as levels) as fixed factor and MouseID (M) as a random factor nested in G.

Accordingly, the lme syntax does not consider S and only consider my fixed variable Genotype (G) and the random variable MouseID (M). Without any Dummy variable, the estimate of the intercept is the mean of the APP group (127.58) and the estimate of G_WT is -112.17 (mean of WT - mean of APP), meaning that the mean of the WT group is 15.41. In other words, the model correctly assign the mean of the APP group to beta0 and the difference between the means to beta1.

The anova run on lme and stored in statsMM2 gives the same F and p vales as the anovan function. In addition, the use of the Satterthwaite approximation is necessary to have the proper df2, namely 10 (indeed, I have a balanced design with five mice for each genotype). Without it, the df2 is equal to the number of observation minues the number of groups of the fixed factors (n. slices-n.genotypes: 24-2=22); it follows that the pVal is much smaller (2.18e-12). In addition, the use of REML is necessary as well to have the same values between the two syntaxes.

If anyone has anything to add and/or point out some fallacy/errors of my reasoning, please post an answer as well