[R-sig-ME] ordinal model and lsmeans mean class computing

[espesser]

Dear All,


I tried to understand how the predicted mean class of an ordinal mixed 
model was computed in lsmeans().
(I want to compute mean classes from the fixed effects only; lsmeans 
included all the random terms)

# from the data wine, I estimated the model

library(ordinal)

  data(wine)

  fm22.clmm= clmm(rating~temp+contact +(1|judge),data=wine, Hess=T)

# get mean classes with lsmeans

lsmeans(fm22.clmm,   ~ temp, mode="mean.class")

  temp mean.class        SE df asymp.LCL asymp.UCL

  cold    2.339608 0.166799 NA   2.012688   2.666528

  warm    3.479936 0.201921 NA   3.084178   3.875694

  

Results are averaged over the levels of: contact

Confidence level used: 0.95

Warning message:

In is.na(x) : is.na() applied to non-(list or vector) of type 'NULL'

I supposed that the predicted mean class is the mean class number (from 
1 to 5) weighted by the estimated probabilities. Therefore I tried:

library(plyr)

wine$fitted = fm22.clmm$fitted

ddply(wine, ~temp, summarise, sum(fitted*as.numeric(rating)) /sum(fitted))

   temp       ..1

1 cold 2.269931

2 warm 3.487144

There are some discrepancies with the results of lsmeans.

I got similar (small) discrepancies for other data and clmm models.

I certainly missed something, and I would really appreciate your help!

Best regards,

Robert


here is the sessionInfo:


sessionInfo()
R version 3.2.2 (2015-08-14)
Platform: x86_64-w64-mingw32/x64 (64-bit)
Running under: Windows 7 x64 (build 7601) Service Pack 1


attached base packages:
[1] stats     graphics  grDevices utils     datasets  methods
[7] base

other attached packages:
[1] plyr_1.8.3         ordinal_2015.6-28  lsmeans_2.23
[4] estimability_1.1-1

loaded via a namespace (and not attached):
  [1] Rcpp_0.12.2      lattice_0.20-33  codetools_0.2-14
  [4] mvtnorm_1.0-5    zoo_1.7-13       ucminf_1.1-3
  [7] MASS_7.3-45      grid_3.2.2       xtable_1.8-2
[10] nlme_3.1-121     coda_0.18-1      multcomp_1.4-5
[13] Matrix_1.2-2     sandwich_2.3-4   splines_3.2.2
[16] TH.data_1.0-7    tools_3.2.2      survival_2.38-3




-- 
Robert Espesser
CNRS UMR  7309 - Université Aix-Marseille
5 Avenue Pasteur
13100 AIX-EN-PROVENCE

Tel: +33 (0)413 55 36 26


	[[alternative HTML version deleted]]