[R] F-test degree of freedoms in lme4 ?

From: Wilhelm B. Kloke <wb_at_arb-phys.uni-dortmund.de>
Date: Thu 12 Jan 2006 - 03:53:25 EST


I have a problem moving from multistratum aov analysis to lmer.

My dataset has observations of ampl at 4 levels of gapf and 2 levels of bl on 6 subjects levels VP, with 2 replicates wg each, and is balanced.

Here is the summary of this set with aov:
>> summary(aov(ampl~gapf*bl+Error(VP/(bl*gapf)),hframe2))
>
>Error: VP
> Df Sum Sq Mean Sq F value Pr(>F)
>Residuals 5 531 106
>
>Error: VP:bl
> Df Sum Sq Mean Sq F value Pr(>F)
>bl 1 1700 1700 37.8 0.0017 **
>Residuals 5 225 45
>---
>Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
>
>Error: VP:gapf
> Df Sum Sq Mean Sq F value Pr(>F)
>gapf 3 933 311 24.2 5.3e-06 ***
>Residuals 15 193 13
>---
>Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
>
>Error: VP:bl:gapf
> Df Sum Sq Mean Sq F value Pr(>F)
>gapf:bl 3 93.9 31.3 3.68 0.036 *
>Residuals 15 127.6 8.5
>---
>Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
>
>Error: Within
> Df Sum Sq Mean Sq F value Pr(>F)
>Residuals 48 318 7
>
This is mostly identical the analysis by BMDP 4V, except for the Greenhouse-Geisser epsilons, which are not estimated this way.

I have to analyse a similar dataset, which is not balanced. So I need to change the method. Following Pinheiro/Bates p.90f, I tried
> hf2.lme <- lme(ampl~gapf*bl,hframe2,random=list(VP=pdDiag(~gapf*bl),bl=pdDiag(~gapf)))
and some variations of this to get the same F tests generated. At least, I got the F-test on error stratum VP:bl this way, but not the other two:
>> anova(hf2.lme)
> numDF denDF F-value p-value
>(Intercept) 1 78 764.86 <.0001
>gapf 3 78 17.68 <.0001
>bl 1 5 37.81 0.0017
>gapf:bl 3 78 2.99 0.0362

Then I tried to move to lmer.
I tried to find something equivalent to the above lme call, with no success at all.

In case, that the problem is in the data, here is the set:

VP ampl wg bl gapf
1 WJ 22 w s 144
2 CR 23 w s 144
3 MZ 25 w s 144
4 MP 34 w s 144
5 HJ 36 w s 144
6 SJ 26 w s 144
7 WJ 34 w s 80
8 CR 31 w s 80
9 MZ 33 w s 80
10 MP 36 w s 80
11 HJ 37 w s 80
12 SJ 32 w s 80
13 WJ 34 w s 48
14 CR 37 w s 48
15 MZ 38 w s 48
16 MP 38 w s 48
17 HJ 40 w s 48
18 SJ 32 w s 48
19 WJ 36 w s 16
20 CR 40 w s 16
21 MZ 39 w s 16
22 MP 40 w s 16
23 HJ 40 w s 16
24 SJ 38 w s 16
25 WJ 16 g s 144
26 CR 28 g s 144
27 MZ 18 g s 144
28 MP 33 g s 144
29 HJ 37 g s 144
30 SJ 28 g s 144
31 WJ 28 g s 80
32 CR 33 g s 80
33 MZ 24 g s 80
34 MP 34 g s 80
35 HJ 36 g s 80
36 SJ 30 g s 80
37 WJ 32 g s 48
38 CR 38 g s 48
39 MZ 34 g s 48
40 MP 37 g s 48
41 HJ 39 g s 48
42 SJ 30 g s 48
43 WJ 36 g s 16
44 CR 34 g s 16
45 MZ 36 g s 16
46 MP 40 g s 16
47 HJ 40 g s 16
48 SJ 36 g s 16
49 WJ 22 w b 144
50 CR 24 w b 144
51 MZ 20 w b 144
52 MP 26 w b 144
53 HJ 22 w b 144
54 SJ 16 w b 144
55 WJ 26 w b 80
56 CR 24 w b 80
57 MZ 26 w b 80
58 MP 27 w b 80
59 HJ 26 w b 80
60 SJ 18 w b 80
61 WJ 28 w b 48
62 CR 23 w b 48
63 MZ 28 w b 48
64 MP 29 w b 48
65 HJ 27 w b 48
66 SJ 24 w b 48
67 WJ 32 w b 16
68 CR 26 w b 16
69 MZ 30 w b 16
70 MP 28 w b 16
71 HJ 30 w b 16
72 SJ 22 w b 16
73 WJ 22 g b 144
74 CR 18 g b 144
75 MZ 18 g b 144
76 MP 26 g b 144
77 HJ 22 g b 144
78 SJ 18 g b 144
79 WJ 24 g b 80
80 CR 26 g b 80
81 MZ 30 g b 80
82 MP 26 g b 80
83 HJ 26 g b 80
84 SJ 24 g b 80
85 WJ 28 g b 48
86 CR 28 g b 48
87 MZ 27 g b 48
88 MP 30 g b 48
89 HJ 26 g b 48
90 SJ 16 g b 48
91 WJ 28 g b 16
92 CR 19 g b 16
93 MZ 24 g b 16
94 MP 32 g b 16
95 HJ 30 g b 16
96 SJ 22 g b 16

-- 
Dipl.-Math. Wilhelm Bernhard Kloke
Institut fuer Arbeitsphysiologie an der Universitaet Dortmund
Ardeystrasse 67, D-44139 Dortmund, Tel. 0231-1084-257

______________________________________________
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Received on Thu Jan 12 08:28:41 2006

This archive was generated by hypermail 2.1.8 : Fri 03 Mar 2006 - 03:42:03 EST