[R-sig-ME] Apparent false convergence in lmer with some BIBD data

[Phillip Chapman]

Greetings ME Board members,

Is it possible to manipulate the convergence parameters in lmer?
The control options on the lmer help page are very limited.  I would like
to alter the convergence criteria.  Also, it would be nice to be able to
input my own starting values. Can someone tell me whether this is 
possible and
point me to an appropriate reference.

I am using the "Pillow" BIBD data from page 1066 of Ott and
Longnecker, 5th edition. My lmer program seems to have a false 
convergence with
the random block estimate at zero.  At the risk of committing heresy, I 
have posted
SAS analysis of the same data.  The SAS solution has a lower 
REMLDeviance than the lmer solution.
I also ran the SAS code with the lmer solution to confirm that the 
problem is with the
optimization, not with the data or the REML likelihood.  lmer and SAS 
give the same
REMLDeviance at the lmer solution.  I have plotted the likelihood 
surface (vertical=block
variance, horizontal=error variance) and see that there is a ridge 
running from SW to NE on
plot.  The lmer solution is not on that ridge.  

Any help or suggestions would be appreciated.

Thanks,
Phil Chapman

Here is the data:

blk pillow firmness
1 A 59
1 B 26
1 C 38
2 D 85
2 E 92
2 F 69
3 G 74
3 H 52
3 I 27
4 A 62
4 D 70
4 G 68
5 B 27
5 E 98
5 H 59
6 C 31
6 F 60
6 I 35
7 A 63
7 E 85
7 I 30
8 B 22
8 F 73
8 G 75
9 C 45
9 D 74
9 H 51
10 A 52
10 F 76
10 H 43
11 B 18
11 D 79
11 I 41
12 C 41
12 E 84
12 G 81

Here is my program and some of its output.

data <- read.table(file="c:/temp/Pillow data.txt",header=T)
attach(data)

blk <- factor(blk)

install.packages("lme4")
library(lme4)

h <- lmer(firmness ~ pillow + (1|blk), REML=TRUE, verbose=TRUE)
print(h, digits=10)

Linear mixed model fit by REML
Formula: firmness ~ pillow + (1 | blk)
         AIC         BIC       logLik    deviance     REMLdev
 207.3991367 224.8178450 -92.69956833 220.2067619 185.3991367
Random effects:
 Groups   Name        Variance   Std.Dev.  
 blk      (Intercept) 2.6584e-10 1.6305e-05
 Residual             3.5398e+01 5.9496e+00
Number of obs: 36, groups: blk, 12




Here is a SAS program with the same model:

proc mixed data=pillow;
class blk pillow;
model firmness=pillow;
random blk;

                           The Mixed Procedure

                            Iteration History

       Iteration    Evaluations    -2 Res Log Like       Criterion

               0              1       185.39913672
               1              2       185.35185521      0.00000016
               2              1       185.35184420      0.00000000


                        Convergence criteria met.


                          Covariance Parameter
                                Estimates

                          Cov Parm     Estimate

                          blk            2.0411
                          Residual      33.3815


                             Fit Statistics

                  -2 Res Log Likelihood           185.4
                  AIC (smaller is better)         189.4
                  AICC (smaller is better)        189.9
                  BIC (smaller is better)         190.3


                      Type 3 Tests of Fixed Effects

                            Num     Den
              Effect         DF      DF    F Value    Pr > F

              pillow          8      16      57.59    <.0001





Here is SAS code with the parameters fixed at the lmer variance estimates:


proc mixed data=pillow;
class blk pillow;
model firmness=pillow;
random blk;
parms (2.6584e-10) (3.5398e+01) /noiter;
run;



                           The Mixed Procedure

                             Parameter Search

    CovP1      CovP2   Variance           Res Log Like   -2 Res Log Like

 2.66E-10    35.3980    35.3981               -92.6996          185.3991


                          Covariance Parameter
                                Estimates

                          Cov Parm     Estimate

                          blk          2.66E-10
                          Residual      35.3981


                             Fit Statistics

                  -2 Res Log Likelihood           185.4
                  AIC (smaller is better)         185.4
                  AICC (smaller is better)        185.4
                  BIC (smaller is better)         185.4


                      Type 3 Tests of Fixed Effects

                            Num     Den
              Effect         DF      DF    F Value    Pr > F

              pillow          8      16      56.60    <.0001