[R] two-level poisson, again

Henric Nilsson henric.nilsson at statisticon.se
Wed Aug 17 10:14:46 CEST 2005


On On, 2005-08-17, 08:52, Shige Song skrev:
> Hi,
>
> I compare results of a simple two-level poisson estimated using lmer
> and those estimated using MLwiN and Stata (v.9).
>
> In R, I trype:
> -------------------------------------------------------------------------------------------
> m2 <- lmer(.D ~ offset(log(.Y)) + (1|pcid2) + educy + agri, male, poisson)
> -------------------------------------------------------------------------------------------
>
> In Stata, I type:
> -------------------------------------------------------------------------------------------
> xtpois _D educy agri, e(_Y) i(pcid2)

You're not fitting the same model! `lmer' uses Gaussian random effects and
the default for `xtpois' is gamma random effects.

Also, note that even if you'd specified a Gaussian random effect (through
a `normal' to the right of the `,' in your `xtpois' call) the same fitting
criterion is not used since `xtpois' uses adaptive Gauss-Hermite
quadrature and `lmer' defaults to PQL.

For comparable results, try the following:

m2 <- lmer(.D ~ offset(log(.Y)) + (1|pcid2) + educy + agri, male, poisson,
method = "AGQ")

xtpois _D educy agri, e(_Y) i(pcid2) re normal

You may also want to try Göran Broström's `glmmML' package.


HTH,
Henric

>
> Results using R look like:
> -------------------------------------------------------------------------------------------
> ..
> Random effects:
>      Groups        Name    Variance    Std.Dev.
>       pcid2 (Intercept)       5e-10  2.2361e-05
> # of obs: 25360, groups: pcid2, 174
>
> Estimated scale (compare to 1)  7.190793
>
> Fixed effects:
>                Estimate  Std. Error z value  Pr(>|z|)
> (Intercept) -3.28548086  0.00408771 -803.75 < 2.2e-16 ***
> educy        0.00507475  0.00039616   12.81 < 2.2e-16 ***
> agri         0.01346887  0.00334837    4.02 5.758e-05 ***
> ..
> ------------------------------------------------------------------------------------------
>
> Results using Stata look like:
>
> ------------------------------------------------------------------------------
>           _D |      Coef.   Std. Err.      z    P>|z|     [95% Conf.
> Interval]
> -------------+----------------------------------------------------------------
>        educy |   .0120431   .0004441    27.12   0.000     .0111725
> .0129136
>         agri |   .0293177   .0035586     8.24   0.000      .022343
> .0362924
>        _cons |  -3.325073   .0076275  -435.93   0.000    -3.340023
> -3.310123
>           _Y | (exposure)
> -------------+----------------------------------------------------------------
>     /lnalpha |  -4.982977   .1156474                     -5.209641
> -4.756312
> -------------+----------------------------------------------------------------
>        alpha |   .0068536   .0007926                      .0054636
> .0085973
> ------------------------------------------------------------------------------
>
>
> As you can see, the discrepency is significant! And results using
> MLwiN agree with Stata. Any help will be greately appreciated!
>
> Shige
>
> ______________________________________________
> R-help at stat.math.ethz.ch mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide!
> http://www.R-project.org/posting-guide.html




More information about the R-help mailing list