Dear lmer-users:
After having reading as many documents as I could about the lmer
function and the lack of p-values, I need to recognize that I haven?t
understood the right way of extracting conclusions when using the lmer
function. Maybe I would need more statistical knowledge or read any
other book (please, any suggestion?), but if you give me any help, I
would be very grateful.
I?m trying to analyze a field experiment with a split-plot design. We
have selected 9 plots of similar area in field: 3 of them were woodland,
3 were shrubland, and the other 3 were open areas. In each plot, we set
up 20 subplots: 10 of them were regularly watered and the other 10 were
control. In each subplot, we sowed seeds of 5 species, each species
protected by a mesh cage. The number of seeds sowed in each cage was
different for each species, ranging from 5 to 15 seeds. So, there are
three fixed factors (habitat, species, and watering) and two random
factors (plot, subplot). I have labelled each plot (n=9) and each
subplot (n=180) with a different number per each variable. I consider
that the seeds sowed in the same cage are not independent, so I used the
counts of emerged/survived seeds per cage as response variable.
Both types of response variables I am interested in (emergence and
surviva) follow a binomial distribution. My question is: is there any
effect of the habitat, the watering treatment and the species identity
for the emergence/survival of seedlings? To answer this question, I
think the best way is using lmer as follows:
emerg1<-lmer(cbind(data$emerged,data$nonemerg)~ habitat * species *
watering + (1|plot) + (1|subplot), data = data, family = binomial)
summary (emerg1)
Generalized linear mixed model fit using Laplace
Formula: cbind(data$emerged, data$nonemerg) ~ habitat * species *
watering + (1 | plot) + (1 | subplot)
Data: data
Family: binomial(logit link)
AIC BIC logLik deviance
2048 2202 -992.2 1984
Random effects:
Groups Name Variance Std.Dev.
subplot (Intercept) 0.08740 0.29564
plot (Intercept) 0.15162 0.38939
number of obs: 900, groups: subplot, 180; plot, 9
Estimated scale (compare to 1 ) 1.341165
Fixed effects:
Estimate Std. Error z value Pr(>|z|)
(Intercept) -2.2724 0.2811 -8.084 6.25e-16 ***
habitatShrubland -0.2837 0.4050 -0.700 0.48363
habitatWoodland -0.1464 0.4005 -0.365 0.71479
speciesencina 2.2170 0.2302 9.630 < 2e-16 ***
speciespino 2.6165 0.1872 13.976 < 2e-16 ***
speciesroble 1.9082 0.2318 8.231 < 2e-16 ***
speciessorbus 0.4327 0.2093 2.068 0.03865 *
wateringW -0.2848 0.2513 -1.133 0.25702
habitatShrubland:speciesencina -0.2105 0.3383 -0.622 0.53371
habitatWoodland:speciesencina 0.1260 0.3302 0.382 0.70272
habitatShrubland:speciespino 0.7979 0.2792 2.858 0.00427 **
habitatWoodland:speciespino -0.7454 0.2698 -2.763 0.00573 **
habitatShrubland:speciesroble 0.2184 0.3383 0.645 0.51868
habitatWoodland:speciesroble 0.1152 0.3331 0.346 0.72954
habitatShrubland:speciessorbus 0.4873 0.3016 1.616 0.10612
habitatWoodland:speciessorbus 0.7417 0.2902 2.556 0.01059 *
habitatShrubland:wateringW 0.4320 0.3588 1.204 0.22866
habitatWoodland:wateringW 0.4770 0.3477 1.372 0.17011
speciesencina:wateringW 0.6484 0.3364 1.927 0.05392 .
speciespino:wateringW 0.2311 0.2762 0.837 0.40277
speciesroble:wateringW 0.1911 0.3383 0.565 0.57220
speciessorbus:wateringW 0.4430 0.3032 1.461 0.14397
habitatShrubland:speciesencina:wateringW -0.5240 0.4806 -1.090 0.27551
habitatWoodland:speciesencina:wateringW -0.6012 0.4703 -1.278 0.20120
habitatShrubland:speciespino:wateringW -0.2279 0.3988 -0.572 0.56764
habitatWoodland:speciespino:wateringW -0.2531 0.3842 -0.659 0.50999
habitatShrubland:speciesroble:wateringW -0.1877 0.4811 -0.390 0.69640
habitatWoodland:speciesroble:wateringW -1.0787 0.4813 -2.241 0.02502 *
habitatShrubland:speciessorbus:wateringW -0.5306 0.4276 -1.241 0.21459
habitatWoodland:speciessorbus:wateringW -0.1205 0.4083 -0.295 0.76791
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
.....
> library(coda)
> A<-mcmcsamp(emerg1,50000)
> summary(A)
Iterations = 1:50000
Thinning interval = 1
Number of chains = 1
Sample size per chain = 50000
1. Empirical mean and standard deviation for each variable,
plus standard error of the mean:
Mean SD Naive SE Time-series SE
(Intercept) -2.2730 0.3903 0.0017456 0.003199
habitatShrubland -0.3034 0.5584 0.0024972 0.003928
habitatWoodland -0.1641 0.5554 0.0024839 0.005055
speciesencina 2.2197 0.2319 0.0010373 0.001756
speciespino 2.6208 0.1879 0.0008401 0.001669
speciesroble 1.9081 0.2323 0.0010390 0.001778
speciessorbus 0.4299 0.2092 0.0009355 0.001784
wateringW -0.2991 0.2504 0.0011199 0.002122
habitatShrubland:speciesencina -0.1947 0.3387 0.0015145 0.002391
habitatWoodland:speciesencina 0.1405 0.3329 0.0014886 0.002777
habitatShrubland:speciespino 0.8186 0.2821 0.0012618 0.002427
habitatWoodland:speciespino -0.7339 0.2710 0.0012118 0.002166
habitatShrubland:speciesroble 0.2349 0.3413 0.0015264 0.002496
habitatWoodland:speciesroble 0.1314 0.3341 0.0014942 0.002882
habitatShrubland:speciessorbus 0.5014 0.3041 0.0013601 0.002638
habitatWoodland:speciessorbus 0.7566 0.2914 0.0013032 0.002556
habitatShrubland:wateringW 0.4585 0.3603 0.0016115 0.003373
habitatWoodland:wateringW 0.5003 0.3461 0.0015480 0.003174
speciesencina:wateringW 0.6698 0.3382 0.0015127 0.002629
speciespino:wateringW 0.2466 0.2746 0.0012283 0.002280
speciesroble:wateringW 0.2063 0.3367 0.0015059 0.002922
speciessorbus:wateringW 0.4587 0.3028 0.0013541 0.002368
habitatShrubland:speciesencina:wateringW -0.5593 0.4800 0.0021466 0.004227
habitatWoodland:speciesencina:wateringW -0.6300 0.4727 0.0021138 0.003775
habitatShrubland:speciespino:wateringW -0.2546 0.3986 0.0017826 0.003616
habitatWoodland:speciespino:wateringW -0.2791 0.3831 0.0017135 0.003225
habitatShrubland:speciesroble:wateringW -0.2116 0.4841 0.0021648 0.004105
habitatWoodland:speciesroble:wateringW -1.1116 0.4825 0.0021577 0.004582
habitatShrubland:speciessorbus:wateringW -0.5540 0.4278 0.0019133 0.003664
habitatWoodland:speciessorbus:wateringW -0.1416 0.4077 0.0018232 0.003570
log(sbpl.(In)) -2.4073 0.2585 0.0011562 0.005921
log(plot.(In)) -1.2932 0.6567 0.0029369 0.006396
2. Quantiles for each variable:
2.5% 25% 50% 75% 97.5%
(Intercept) -3.04171 -2.506810 -2.2708 -2.03446 -1.5193
habitatShrubland -1.38115 -0.646627 -0.3015 0.02902 0.7991
habitatWoodland -1.26227 -0.501688 -0.1663 0.16544 0.9489
speciesencina 1.77105 2.060916 2.2197 2.37579 2.6798
speciespino 2.26105 2.495121 2.6158 2.74363 3.0019
speciesroble 1.45680 1.751212 1.9071 2.06305 2.3624
speciessorbus 0.02693 0.287434 0.4281 0.56880 0.8456
wateringW -0.79641 -0.469682 -0.2973 -0.12805 0.1822
habitatShrubland:speciesencina -0.85631 -0.424483 -0.1910 0.03690 0.4640
habitatWoodland:speciesencina -0.51369 -0.081558 0.1417 0.36531 0.7943
habitatShrubland:speciespino 0.26948 0.628067 0.8161 1.00870 1.3755
habitatWoodland:speciespino -1.26446 -0.916398 -0.7333 -0.55441 -0.2018
habitatShrubland:speciesroble -0.43404 0.002914 0.2329 0.46687 0.9052
habitatWoodland:speciesroble -0.51401 -0.098191 0.1332 0.35979 0.7730
habitatShrubland:speciessorbus -0.09194 0.295077 0.5010 0.70509 1.0973
habitatWoodland:speciessorbus 0.18432 0.558907 0.7553 0.95078 1.3234
habitatShrubland:wateringW -0.25126 0.218858 0.4567 0.69857 1.1727
habitatWoodland:wateringW -0.17991 0.264459 0.5001 0.72907 1.1924
speciesencina:wateringW 0.01156 0.437409 0.6705 0.89904 1.3409
speciespino:wateringW -0.28788 0.059775 0.2445 0.43289 0.7807
speciesroble:wateringW -0.44262 -0.021072 0.2011 0.42942 0.8769
speciessorbus:wateringW -0.12375 0.248090 0.4556 0.66512 1.0558
habitatShrubland:speciesencina:wateringW -1.50154 -0.886472 -0.5554
-0.23413 0.3715
habitatWoodland:speciesencina:wateringW -1.55928 -0.942517 -0.6263
-0.31732 0.2906
habitatShrubland:speciespino:wateringW -1.02972 -0.524534 -0.2583
0.01483 0.5318
habitatWoodland:speciespino:wateringW -1.03432 -0.536075 -0.2782
-0.02089 0.4624
habitatShrubland:speciesroble:wateringW -1.16509 -0.538283 -0.2073
0.10929 0.7355
habitatWoodland:speciesroble:wateringW -2.05038 -1.439554 -1.1132
-0.78493 -0.1720
habitatShrubland:speciessorbus:wateringW -1.39134 -0.841315 -0.5548
-0.26179 0.2815
habitatWoodland:speciessorbus:wateringW -0.94322 -0.415447 -0.1360
0.13537 0.6568
log(sbpl.(In)) -2.95653 -2.567736 -2.3935 -2.22955 -1.9423
log(plot.(In)) -2.41148 -1.753486 -1.3534 -0.89638 0.1613
> HPDinterval(A)
lower upper
(Intercept) -3.02693912 -1.50685625
habitatShrubland -1.37402561 0.80185848
habitatWoodland -1.28262371 0.92217838
speciesencina 1.77588054 2.68156075
speciespino 2.24528494 2.98347131
speciesroble 1.45458293 2.35890008
speciessorbus 0.02277444 0.83992583
wateringW -0.79840382 0.17786065
habitatShrubland:speciesencina -0.87303659 0.44080390
habitatWoodland:speciesencina -0.47325939 0.82869813
habitatShrubland:speciespino 0.26283252 1.36579201
habitatWoodland:speciespino -1.26969655 -0.20969532
habitatShrubland:speciesroble -0.44787916 0.89012541
habitatWoodland:speciesroble -0.50423916 0.78148583
habitatShrubland:speciessorbus -0.09267742 1.09603954
habitatWoodland:speciessorbus 0.18433471 1.32351169
habitatShrubland:wateringW -0.26481667 1.15194682
habitatWoodland:wateringW -0.15799665 1.20945912
speciesencina:wateringW 0.02501399 1.34820459
speciespino:wateringW -0.29087519 0.77543688
speciesroble:wateringW -0.44362264 0.87502474
speciessorbus:wateringW -0.11402239 1.06217002
habitatShrubland:speciesencina:wateringW -1.48872380 0.38175569
habitatWoodland:speciesencina:wateringW -1.54261721 0.30557576
habitatShrubland:speciespino:wateringW -1.01824049 0.53687967
habitatWoodland:speciespino:wateringW -1.03432105 0.46236822
habitatShrubland:speciesroble:wateringW -1.17560521 0.72330469
habitatWoodland:speciesroble:wateringW -2.05138457 -0.17460151
habitatShrubland:speciessorbus:wateringW -1.39314261 0.27800441
habitatWoodland:speciessorbus:wateringW -0.93252386 0.66400031
log(sbpl.(In)) -2.92059137 -1.91602692
log(plot.(In)) -2.49817299 0.03793472
attr(,"Probability")
[1] 0.95
>
Sincerely, for me now is the great deal. Without any p-value, how can I
know if a fixed factor is significant for the response variable? What is
the interpretation of these results? How I should present the
information in a correct way to editors and referees (and for me, to
understanding the effects)?
I am afraid you maybe have answered this a thousand times, but I have
read the wiki about lmer
(http://wiki.r-project.org/rwiki/doku.php?id=guides:lmer-tests&s=lme%20and%20aov)
and the discussion about the p-values, and I still feel confused with
the interpretation of results. Do you still think that it is correct the
mcmcpvalue function written by Douglas Bates? Any help is welcome.
Thanks a lot!
Greetings,
Irene
Irene Mendoza Sagrera
PhD Student
Terrestrial Ecology Group.
University of Granada (Spain)
Departamento de Ecolog?a
Facultad de Ciencias.
Campus de Fuentenueva, s/n.
Universidad de Granada, 18071.
Granada, Spain
Phone: +34 958 243242
Fax: +34 958 243238
website www.ugr.es/~rnm220
website www.ugr.es/~redbome
After having reading as many documents as I could about the lmer
function and the lack of p-values, I need to recognize that I
haven?t understood the right way of extracting conclusions when
using the lmer function. Maybe I would need more statistical
knowledge or read any other book (please, any suggestion?), but if
you give me any help, I would be very grateful.
I?m trying to analyze a field experiment with a split-plot design.
...
Sincerely, for me now is the great deal. Without any p-value, how
can I know if a fixed factor is significant for the response
variable? What is the interpretation of these results? How I should
present the information in a correct way to editors and referees
(and for me, to understanding the effects)?
Irene--
While the examples won't be directly related to your field, I have
found that a paper by Baayen, Davidson, and Bates has a very clear
presentation of the analysis of a split-plot design with lme4. In
addition, the article has a corresponding R package called languageR
that provides some functionality that would be useful in answering
your main question. Taken together, I think the paper and the R
package are a very useful addition to the available materials on mixed
effects models. I'm grateful to the authors for putting both
together!
You can find the submitted pdf on Baayen's publications page, at
http://www.mpi.nl/world/persons/private/baayen/publications/baayenDavidsonBates.pdf.
Good luck with your analysis!
/au