nlme estima una varianza cercana a cero para los efectos aleatorios

Question

Estoy realizando diversos análisis sobre una pequeña muestra. Básicamente, tenemos un experimento en el que 14 sujetos (UID 1 ~ 14) utilizaron uno de los 6 instrumentos (MID 1 ~ 6) en 3 ocasiones (Secuencia 1 ~ 3). Cada vez se registró una puntuación de resultado (entre 1 ~ 100).

La prueba fue a doble ciego. Se dijo a los sujetos que estaban midiendo 3 condiciones diferentes, mientras que en realidad estaban midiendo las condiciones A, B, A o B, A, B (asignadas aleatoriamente a las máquinas y a los usuarios). El objetivo era comprobar si A y B son diferentes o no.

Para ver si hay alguna diferencia significativa entre las puntuaciones de las condiciones A y B, he intentado ajustar un modelo de intercepción aleatoria simple utilizando el paquete nlme de R. He probado:

f.1 <- lme(Score ~ Condition, random = ~1|UID, data)

Sin embargo, por alguna razón lme no se ajusta al modelo: no da ningún error ni advertencia, pero la varianza del efecto aleatorio ajustado es esencialmente cero:

> summary(f.1)
Linear mixed-effects model fit by REML
 Data: data 
       AIC      BIC   logLik
  349.3259 356.0815 -170.663

Random effects:
 Formula: ~1 | UID
         (Intercept) Residual
StdDev: 0.0009303203 15.98295

Fixed effects: Score ~ Condition 
               Value Std.Error DF   t-value p-value
(Intercept) 77.47619  3.487766 27 22.213700  0.0000
ConditionA  -0.85714  4.932446 27 -0.173776  0.8633
 Correlation: 
           (Intr)
ConditionA -0.707

Standardized Within-Group Residuals:
       Min         Q1        Med         Q3        Max 
-2.9704269 -0.4677603  0.2472873  0.7835730  1.4628682 

Number of Observations: 42
Number of Groups: 14

He intentado hacer lo mismo utilizando lme4 y obtuve los mismos resultados. Las estimaciones para el intercepto y el Condition factor es casi idéntico a un modelo lineal si utilizo lm .

Me esfuerzo por comprender lo que lme o lmer fallan al estimar el efecto aleatorio. He generado algunos datos por simulación y ambas rutinas no han tenido problemas para ajustar el modelo, por lo que dudo que haya algo mal en la sintaxis de lo que he utilizado.

Los datos están aquí:

   UID MID Seq Score Condition
1    1   1   1    90  B
2    1   1   2    85  A
3    1   1   3    75  B
4    2   4   1    75  A
5    2   4   2    95  B
6    2   4   3    85  A
7    3   6   1    60  A
8    3   6   2    82  B
9    3   6   3    85  A
10   4   3   1    60  A
11   4   3   2    70  B
12   4   3   3    75  A
13   5   2   1    85  B
14   5   2   2    85  A
15   5   2   3    85  B
16   6   5   1    90  B
17   6   5   2    95  A
18   6   5   3   100  B
19   7   2   1    90  B
20   7   2   2    70  A
21   7   2   3    50  B
22   8   1   1    70  B
23   8   1   2    75  A
24   8   1   3    80  B
25   9   3   1    90  A
26   9   3   2    30  B
27   9   3   3    90  A
28  10   6   1    50  A
29  10   6   2    85  B
30  10   6   3    92  A
31  11   4   1    50  A
32  11   4   2    85  B
33  11   4   3    92  A
34  12   5   1    65  B
35  12   5   2    50  A
36  12   5   3    90  B
37  13   4   1    65  A
38  13   4   2    70  B
39  13   4   3    80  A
40  14   2   1    60  B
41  14   2   2   100  A
42  14   2   3    80  B

Answer 1

Como se menciona en los comentarios a la pregunta, el problema es un ajuste singuar. Una solución es adoptar un enfoque bayesiano, que aparentemente ha sido una buena solución en este caso.