Hi All

Just to follow up on my original posting. I checked the results against lme( ) (from the nlme library) and now I am getting identical results to asreml( ) whereas lmer( ) from the lme4 library gives the wrong values for the standard errors and estimates of the variance components (see below). Note that I have now corrected an earlier fault where I failed to correctly scale standard errors by sigma in the asreml output. So there is definitely a fault with lmer( ) and users should be aware of this problem.

> lme.01 <- lmer(log10Abund ~ cell.f + YEAR.vcentred + (1|YEAR.fs), weights=data.tab.sci$vT_noshots, data=data.tab.sci)
>
> summary(lme.01)
Linear mixed model fit by REML
Formula: log10Abund ~ cell.f + YEAR.vcentred + (1 | YEAR.fs)
   Data: data.tab.sci
   AIC   BIC logLik deviance REMLdev
542.7 577.3 -258.4    496.8   516.7
Random effects:
Groups   Name        Variance  Std.Dev.
YEAR.fs  (Intercept) 0.0041281 0.06425
 Residual             0.1630013 0.40373
Number of obs: 106, groups: YEAR.fs, 23

Fixed effects:
               Estimate Std. Error t value
(Intercept)    0.587934   0.046173  12.733
cell.f202      0.471476   0.053572   8.801
cell.f203      0.563821   0.045650  12.351
cell.f303      0.408994   0.051588   7.928
cell.f304      0.004537   0.066063   0.069
cell.f403      0.329375   0.076186   4.323
cell.f404      0.390589   0.072599   5.380
cell.f405     -0.032795   0.066021  -0.497
cell.f406     -0.153596   0.054923  -2.797
cell.f505      0.979260   0.059209  16.539
YEAR.vcentred -0.013434   0.002012  -6.676

[snip]
> data.tab.sci$IvT_noshots <- 1/data.tab.sci$vT_noshots

> nlme.01 <- lme(fixed=log10Abund ~ cell.f + YEAR.vcentred, random= ~1|YEAR.fs, weights=varFixed(~IvT_noshots), data=data.tab.sci)
>
> summary(nlme.01)
Linear mixed-effects model fit by REML
Data: data.tab.sci
       AIC     BIC    logLik
  266.1846 299.385 -120.0923

Random effects:
Formula: ~1 | YEAR.fs
        (Intercept) Residual
StdDev:    0.363026 2.281171

Variance function:
Structure: fixed weights
Formula: ~IvT_noshots
Fixed effects: log10Abund ~ cell.f + YEAR.vcentred
                   Value Std.Error DF    t-value p-value
(Intercept)    0.5879281 0.2608859 74  2.2535835  0.0272
cell.f202      0.4714769 0.3026892 74  1.5576270  0.1236
cell.f203      0.5638270 0.2579314 74  2.1859567  0.0320
cell.f303      0.4089995 0.2914852 74  1.4031570  0.1648
cell.f304      0.0045382 0.3732673 74  0.0121580  0.9903
cell.f403      0.3293774 0.4304649 74  0.7651667  0.4466
cell.f404      0.3905906 0.4101973 74  0.9522018  0.3441
cell.f405     -0.0327946 0.3730345 74 -0.0879131  0.9302
cell.f406     -0.1535953 0.3103284 74 -0.4949443  0.6221
cell.f505      0.9792600 0.3345431 74  2.9271567  0.0045
YEAR.vcentred -0.0134340 0.0113705 21 -1.1814839  0.2506
[snip]

> asreml.01 <-asreml(fixed = log10Abund ~ cell.f + YEAR.vcentred,
+        random=~YEAR.fs, weights=vT_noshots,
+        data=data.tab.sci, na.method.Y = "exclude",
+        na.method.X = "exclude", maxiter=50)

asreml(): 3.0.1 Library: 3.01gl IA32  Run: Mon Feb 25 12:00:24 2013

     LogLik         S2      DF
    -36.9734      4.4391    95  12:00:24
    -33.0962      4.9321    95  12:00:24
    -32.8240      5.1107    95  12:00:24
    -32.7934      5.2124    95  12:00:24
    -32.7932      5.2037    95  12:00:24
    -32.7932      5.2037    95  12:00:24

Finished on: Mon Feb 25 12:00:24 2013
LogLikelihood Converged
>
> summary(asreml.01)
$call
asreml(fixed = log10Abund ~ cell.f + YEAR.vcentred, random = ~YEAR.fs,
    data = data.tab.sci, weights = vT_noshots, na.method.Y = "exclude",
    na.method.X = "exclude", maxiter = 50)

$loglik
[1] -32.79315

$nedf
[1] 95

$sigma
[1] 2.281171

$varcomp
                         gamma component  std.error  z.ratio constraint
YEAR.fs!YEAR.fs.var 0.02532561 0.1317879 0.05610776 2.348835   Positive
R!variance          1.00000000 5.2037401 0.84908182 6.128667   Positive

attr(,"class")
[1] "summary.asreml"
>
> asreml.01$coefficients$fixed
YEAR.vcentred    cell.f_102    cell.f_202    cell.f_203    cell.f_303
 -0.013434004   0.000000000   0.471476899   0.563826976   0.408999494
   cell.f_304    cell.f_403    cell.f_404    cell.f_405    cell.f_406
  0.004538179   0.329377419   0.390590617  -0.032794632  -0.153595263
   cell.f_505   (Intercept)
  0.979259954   0.587928074
> (summary(asreml.01))$sigma*((asreml.01$vcoeff$fixed)^0.5)
[1] 0.01137045 0.00000000 0.30268922 0.25793145 0.29148521 0.37326733
[7] 0.43046490 0.41019732 0.37303451 0.31032841 0.33454305 0.26088586
> asreml.01$coefficients$fixed["YEAR.vcentred"]
YEAR.vcentred
    -0.013434
> (summary(asreml.01))$sigma*((asreml.01$vcoeff$fixed)^0.5)[1]
[1] 0.01137045
>
>
> summary(asreml.01)$varcomp
                         gamma component  std.error  z.ratio constraint
YEAR.fs!YEAR.fs.var 0.02532561 0.1317879 0.05610776 2.348835   Positive
R!variance          1.00000000 5.2037401 0.84908182 6.128667   Positive
> (summary(asreml.01)$varcomp[["component"]][1])^0.5
[1] 0.363026

Information on package 'nlme'

Description:

Package:            nlme
Version:            3.1-104
Date:               2012-05-21

Information on package 'lme4'

Description:

Package:            lme4
Version:            0.999999-0
Date:               2012-06-22



___________________________________________________________________________

    Australian Antarctic Division - Commonwealth of Australia
IMPORTANT: This transmission is intended for the addressee only. If you are not the
intended recipient, you are notified that use or dissemination of this communication is
strictly prohibited by Commonwealth law. If you have received this transmission in error,
please notify the sender immediately by e-mail or by telephoning +61 3 6232 3209 and
DELETE the message.
        Visit our web site at http://www.antarctica.gov.au/
___________________________________________________________________________

	[[alternative HTML version deleted]]

