Mixed effect model with no grouping/cluster

22 views
Skip to first unread message

Pedro Emmanuel Alvarenga Americano do Brasil

unread,
May 20, 2013, 11:50:16 AM5/20/13
to MedStats-list
StatMasters,

One colleague of mine and I are trying to reproduce an example (and we are stuck) from the following paper:

Haobo Ren , Xiao-Hua Zhou & Hua Liang (2004): A Flexible Method for Estimating the ROC Curve, Journal of Applied Statistics, 31:7, 773-784

Here the author argues that it is possible to consistently estimate a robust ROC curve using a mixed effect model. At a certain point the author states that the procedure is easy with a S-PLUS function for mixed effect model - although we are using R, whit the lme function. It seems that the author used the random effect without stating a cluster or grouping information. However, we are kind of banging our heads around for a while as it does not make much sense for us not having a "gouping" variable. 

The author states that the model expression is something similar to:

Rempi = X*Beta + Z*b + e 

where in our example ...

Rempi <- c(0.75, 0.80, 0.85, 0.85, 0.85, 0.85, 0.85, 0.85, 0.85, 0.85, 0.85, 0.90, 0.95, 0.95, 0.95, 1.00, 1.00, 1.00, 1.00)

X <- matrix(c(1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 
1, 1, 0.05, 0.1, 0.15, 0.2, 0.25, 0.3, 0.35, 0.4, 0.45, 0.5, 
0.55, 0.6, 0.65, 0.7, 0.75, 0.8, 0.85, 0.9, 0.95, 0.0025, 0.01, 
0.0225, 0.04, 0.0625, 0.09, 0.1225, 0.16, 0.2025, 0.25, 0.3025, 
0.36, 0.4225, 0.49, 0.5625, 0.64, 0.7225, 0.81, 0.9025, 0.000125, 
0.001, 0.003375, 0.008, 0.015625, 0.027, 0.042875, 0.064, 0.091125, 
0.125, 0.166375, 0.216, 0.274625, 0.343, 0.421875, 0.512, 0.614125, 
0.729, 0.857375),nrow=19,byrow=F)

Z <- matrix(c(0, 0, 2e-04, 0.0012, 0.0037, 0.0086, 0.0166, 0.0284, 
0.0447, 0.0664, 0.0942, 0.1288, 0.171, 0.2214, 0.281, 0.3504, 
0.4304, 0.5217, 0.625, 0, 0, 0, 0, 3e-04, 0.0015, 0.0045, 0.0099, 
0.0186, 0.0313, 0.0486, 0.0715, 0.1005, 0.1366, 0.1804, 0.2326, 
0.2941, 0.3655, 0.4477, 0, 0, 0, 0, 0, 0, 4e-04, 0.002, 0.0054, 
0.0114, 0.0208, 0.0343, 0.0527, 0.0768, 0.1072, 0.1447, 0.1901, 
0.2441, 0.3075, 0, 0, 0, 0, 0, 0, 0, 0, 6e-04, 0.0025, 0.0063, 
0.013, 0.0231, 0.0376, 0.0571, 0.0823, 0.1141, 0.1531, 0.2002, 
0, 0, 0, 0, 0, 0, 0, 0, 0, 1e-04, 9e-04, 0.003, 0.0074, 0.0147, 
0.0257, 0.0411, 0.0616, 0.0881, 0.1213, 0, 0, 0, 0, 0, 0, 0, 
0, 0, 0, 0, 2e-04, 0.0012, 0.0037, 0.0086, 0.0166, 0.0284, 0.0447, 
0.0664, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 3e-04, 0.0015, 
0.0045, 0.0099, 0.0186, 0.0313, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 
0, 0, 0, 0, 0, 4e-04, 0.002, 0.0054, 0.0114, 0, 0, 0, 0, 0, 0, 
0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 6e-04, 0.0025, 0, 0, 0, 0, 0, 
0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1e-04), nrow=19,byrow=F)

In one of our attempts we tried to make each line a single group, but of course it doesnt work. 

> grouping <- 1:19
> lme(fixed = Rempi ~ X, random = ~ Z|grouping)

Erro em MEEM(object, conLin, control$niterEM) : 
  Singularity in backsolve at level 0, block 1
Além disso: Mensagens de aviso perdidas:
In lme.formula(fixed = Rempi ~ X, random = ~Z | grouping) :
  Fewer observations than random effects in all level 1 groups 

Of course one of the things that we tried was to ask the author the original script but no answer came so far. 

We would like to make it work, but the question is: does it make sense to run such analysis without a grouping variable? Any tip to make it work?

Kind regards to all,   

Dr. Pedro Emmanuel A. A. do Brasil
Instituto de Pesquisa Clínica Evandro Chagas
Fundação Oswaldo Cruz
Rio de Janeiro - Brasil
Av. Brasil 4365,
CEP 21040-360,
Tel 55 21 3865-9648
email: pedro....@ipec.fiocruz.br
email: emmanue...@gmail.com

---Apoio aos softwares livres
www.zotero.org - gerenciamento de referências bibliográficas.
www.broffice.org ou www.libreoffice.org - textos, planilhas ou apresentações.
www.epidata.dk - entrada de dados.
www.r-project.org - análise de dados.
www.ubuntu.com - sistema operacional
Reply all
Reply to author
Forward
0 new messages