[R-sig-ME] interpretation of split-plot design with results from lmer

[Irene Mendoza Sagrera]

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.



>  str(data)

'data.frame': 900 obs. of 7 variables:

$ plot : Factor w/ 9 levels "B1","B2","B3",..: 1 1 1 1 1 1 1 1 1 1 ...

$ habitat : Factor w/ 3 levels "Open","Shrubland",..: 1 1 1 1 1 1 1 1 1 
1 ...

$ watering : Factor w/ 2 levels "C","W": 2 2 2 2 2 2 2 2 2 2 ...

$ subplot : Factor w/ 180 levels "1","2","3","4",..: 1 1 1 1 1 2 2 2 2 2 ...

$ species : Factor w/ 5 levels "arce","encina",..: 1 2 3 4 5 1 2 3 4 5 ...

$ emerged : int 3 2 6 1 3 0 3 1 2 2 ...

$ nonemerged: int 12 3 9 4 12 15 2 14 3 13 ...

>



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