[R-sig-ME] naive question about output from mcmcsamp v. lmer estimates

[Gibbons,James]

Dear List,

I am  attempting to estimate the quantitative trait statistic Qst for a 
colleague. Essentially I am fitting the model

X[ijk] = mu + alpha[i] + beta[ij] + gamma[k] + epsilon[ijk]

to data from measured from trees planted in an unbalanced, incomplete 
block experiment, where X[ijk] = individual observation in the ith 
provenance, ijth family and kth block; mu, alpha[i], beta[ij], and 
gamma[k] are the effects of population mean, provenance, family within 
provenance and block respectively. Qst is estimated using variance 
components for provenance and family within provenance. Because of the 
problems inherent in using point estimates to estimate Qst I would like 
to use a Bayesian approach as proposed by Waldmann et al. (Heredity 
(2005) 94, 623–629). I have implemented this model in winBugs and get 
similar summary output to fitting the following lmer model (from lme4, 
using R 2.4.1 on windows):

 >summary(lme.HT02)

Linear mixed-effects model fit by REML
Formula: HT02 ~ (1 | Block) + (1 | Provenance/Family)
    Data: quant
    AIC   BIC logLik MLdeviance REMLdeviance
  12615 12636  -6303      12611        12607
Random effects:
  Groups            Name        Variance Std.Dev.
  Family:Provenance (Intercept)   3.0350  1.7421
  Provenance        (Intercept) 191.7605 13.8478
  Block             (Intercept)  16.1534  4.0191
  Residual                      593.3983 24.3598
number of obs: 1359, groups: Family:Provenance, 335; Provenance, 18; 
Block, 14

Fixed effects:
             Estimate Std. Error t value
(Intercept)  123.497      3.549    34.8

As winBugs is very slow for this problem I was hoping to use mcmcsamp 
instead. I am being naive in thinking that I can use mcmcsamp to 
estimate posterior variance components densities? E.g.

 >colMeans(mcmcsamp(lme.HT02,10000,trans=FALSE))
(Intercept)     sigma^2   Fm:P.(In)   Prvn.(In)   Blck.(In)
   122.57740   571.59628    60.77002    61.59250    17.56870

I was expecting these numbers to be similar, although not identical, to 
the lmer output. This is all quite new ground for me so, I readily 
accept I have probably gone wrong somewhere. But where?

Thanks,

James

-- 
Dr James Gibbons
Research Lecturer in Ecological Modelling
School of the Environment & Natural Resources
University of Wales, Bangor
phone: 01248 382461
email: j.gibbons at bangor.ac.uk


-- 
Gall y neges e-bost hon, ac unrhyw atodiadau a anfonwyd gyda hi,
gynnwys deunydd cyfrinachol ac wedi eu bwriadu i'w defnyddio'n unig
gan y sawl y cawsant eu cyfeirio ato (atynt). Os ydych wedi derbyn y
neges e-bost hon trwy gamgymeriad, rhowch wybod i'r anfonwr ar
unwaith a dilëwch y neges. Os na fwriadwyd anfon y neges atoch chi,
rhaid i chi beidio â defnyddio, cadw neu ddatgelu unrhyw wybodaeth a
gynhwysir ynddi. Mae unrhyw farn neu safbwynt yn eiddo i'r sawl a'i
hanfonodd yn unig  ac nid yw o anghenraid yn cynrychioli barn
Prifysgol Cymru, Bangor. Nid yw Prifysgol Cymru, Bangor yn gwarantu
bod y neges e-bost hon neu unrhyw atodiadau yn rhydd rhag firysau neu
100% yn ddiogel. Oni bai fod hyn wedi ei ddatgan yn uniongyrchol yn
nhestun yr e-bost, nid bwriad y neges e-bost hon yw ffurfio contract
rhwymol - mae rhestr o lofnodwyr awdurdodedig ar gael o Swyddfa
Cyllid Prifysgol Cymru, Bangor.  www.bangor.ac.uk

This email and any attachments may contain confidential mate...{{dropped}}