[R] behaviour of logLik and lme

Dr. Peter Schlattmann peter.schlattmann at charite.de
Tue Apr 19 12:16:33 CEST 2005


Dear all,


when performing a meta analysis I have two results obtained with logLik
and lme, which I do not quite understand.

The results are based on these data:

study  or      var
1   0.10436 0.299111
2  -0.03046 0.121392
3   0.76547 0.319547
4  -0.19845 0.025400
5  -0.10536 0.025041
6  -0.11653 0.040469
7   0.09531 0.026399
8   0.26236 0.017918
9  -0.26136 0.020901
10  0.45742 0.035877
11 -0.59784 0.076356
12 -0.35667 0.186879
13 -0.10536 0.089935
14 -0.31471 0.013772
15 -0.10536 0.089935
16  0.02956 0.004738
17  0.60977 0.035781
18 -0.30111 0.036069
19  0.01980 0.024611
20  0.00000 0.002890
21 -0.04082 0.015863
22  0.02956 0.067069
23  0.18232 0.010677
24  0.26236 0.017918
25  0.32208 0.073896
26  0.67803 0.489415
27 -0.96758 0.194768
28  0.91629 0.051846
29  0.32208 0.110179
30 -1.13943 0.086173
31 -0.47804 0.103522
32  0.16551 0.004152
33  0.46373 0.023150
34 -0.52763 0.050384
35  0.10436 0.003407
36  0.55389 0.054740


The first result concerns logLik

m0<-glm(or~1,family=gaussian(),data=temp,weights=1./var)

logLik(m0)
`log Lik.' -7.10697 (df=2)

For comparison direct calculation of the log likelihood gives:

ll<-sum(log(dnorm(temp$or,fitted(m0),sqrt(temp$var))))
> ll
[1] -33.19137

Does logLik omit constants or how can this discrepancy be explained?


 My second problem is with lme. I want to use FIXED variances in the
estimation process:

m.lm1<-lme(or~1,random=~1|study,weights=varFixed(~var),data=temp,method="ML")


> m.lm1
Linear mixed-effects model fit by maximum likelihood
  Data: temp
  Log-likelihood: -14.22718
  Fixed: or ~ 1
(Intercept)
 0.05597874

Random effects:
 Formula: ~1 | study
        (Intercept) Residual
StdDev: 0.002656846 1.795564

Variance function:
 Structure: fixed weights
 Formula: ~var
Number of Observations: 36
Number of Groups: 36


This is very much the fixed effects result and quite different from the
result I obtain with SAS. Are there any tips to correct this result?

Many thanks in advance

Peter




More information about the R-help mailing list