Sorry,
Actually I gave my data in an image file (.RData) - I've just checked my
send emails.
Am I to give data in another format, such as text ? Here are they in text
(.txt).
The output are :
> summary(aov1 <- aov(Rendement ~ Error(Pollinisateur * Lignee), data =
mca2)
Error: Pollinisateur
Df Sum Sq Mean Sq F value Pr(>F)
Residuals 9 11.9729 1.3303
Error: Lignee
Df Sum Sq Mean Sq F value Pr(>F)
Residuals 4 18.0294 4.5074
Error: Pollinisateur:Lignee
Df Sum Sq Mean Sq F value Pr(>F)
Residuals 36 5.1726 0.1437
Error: Within
Df Sum Sq Mean Sq F value Pr(>F)
Residuals 50 3.7950 0.0759
# F tests :
> Femp <- c(tab1[1:3, 3]/tab1[c(3,3,4), 3])
> names(Femp) <- c("Pollinisateur", "Lignee",
"Interaction")
> Femp
Pollinisateur Lignee Interaction
9.258709 31.370027 1.893061
> 1 - pf(Femp, tab1[1:3,1], tab1[c(3,3,4),1])
Pollinisateur Lignee Interaction
4.230265e-07 2.773448e-11 1.841028e-02
# Standard deviation :
> variances <- c(c(tab1[1:3, 3] - tab1[c(3,3,4), 3]) / c(2*5, 2*10, 2),
tab1[4,3])
> names(variances) <- c(names(Femp), "Residuelle")
> variances
Pollinisateur Lignee Interaction Residuelle
0.11866389 0.21818333 0.03389167 0.07590000
# Using lmer :
> library(lme4)
> summary(lmer(Rendement ~ (1 |Pollinisateur) + (1 | Lignee) + (1 |
Pollinisateur:Lignee),
data=mca2))
Linear mixed-effects model fit by REML
Formula: Rendement ~ (1 | Pollinisateur) + (1 | Lignee) + (1 |
Pollinisateur:Lignee)
Data: mca2
AIC BIC logLik MLdeviance REMLdeviance
105.3845 118.4104 -47.69227 94.35162 95.38453
Random effects:
Groups Name Variance Std.Dev.
Pollinisateur:Lignee (Intercept) 0.033892 0.18410
Pollinisateur (Intercept) 0.118664 0.34448
Lignee (Intercept) 0.218183 0.46710
Residual 0.075900 0.27550
# of obs: 100, groups: Pollinisateur:Lignee, 50; Pollinisateur, 10; Lignee, 5
Fixed effects:
Estimate Std. Error DF t value Pr(>|t|)
(Intercept) 12.60100 0.23862 99 52.808 < 2.2e-16 ***
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05
'.' 0.1 ' ' 1
Thanks,
Jacques VESLOT
Prof Brian Ripley a ?crit :> Nothing was enclosed, nor was the output from summary.aov, so we are
> left guessing.
>
> On Thu, 27 Oct 2005, Jacques VESLOT wrote:
>
>> Dear R-users,
>>
>> My question is how to get right F tests for random effects in random
>> effect models (I hope this question has not been answered too many
>> times yet - I didn't find an answer in rhelp archives).
>>
>> My data are in mca2 (enc.) :
>>
>> names(mca2)
>> [1] "Lignee" "Pollinisateur"
"Rendement"
>>
>> dim(mca2)
>> [1] 100 3
>>
>> replications(Rendement ~ Lignee * Pollinisateur, data = mca2)
>> Lignee Pollinisateur Lignee:Pollinisateur
>> 20 10 2
>>
>> Of course, summary(aov(Rendement ~ Pollinisateur * Lignee, data =
>> mca2)) gives wrong tests of random effects. But, summary(aov1 <-
>> aov(Rendement ~ Error(Pollinisateur * Lignee), data = mca2)) gives no
>> test at all, and I have to do it like this :
>>
>> tab1 <- matrix(unlist(summary(aov1)), nc=5, byrow=T)[,1:3]
>>
>> Femp <- c(tab1[1:3, 3]/tab1[c(3,3,4), 3])
>>
>> names(Femp) <- c("Pollinisateur", "Lignee",
"Interaction")
>>
>> 1 - pf(Femp, tab1[1:3,1], tab1[c(3,3,4),1])
>>
>> With "lme4" package (I did'nt succeed in writing a
working formula
>> with lme from "nlme" package), I can "see" standard
deviations of
>> random effects (but don't know how to find them) with :
>>
>> library(lme4)
>> summary(lmer(Rendement ~ (1 |Pollinisateur) + (1 | Lignee) + (1 |
>> Pollinisateur:Lignee), data=mca2))
>>
>> but I can't get F tests.
>>
>> Thanks in advance.
>>
>> Best regards,
>>
>> Jacques VESLOT
>>
>>
>>
>>
>
-------------- next part --------------
An embedded and charset-unspecified text was scrubbed...
Name: mca2.txt
Url:
https://stat.ethz.ch/pipermail/r-help/attachments/20051027/ce7d42e4/mca2.txt