[ADMB Users] example

[dave fournier]

I see that at the ADMB meeting there was discussion about examples etc.
Here is an example from the R list. It is supposed to show that a
linear mixed model
with only 4 levels of the random effect is a bad idea, ostensibly
because lmer
does not perform well. You can find the discussion at.

https://stat.ethz.ch/pipermail/r-sig-mixed-models/2010q2/003713.html

I think the main point is at
Random effects:
Groups Name Variance Std.Dev.
fac (Intercept) 1.1162 1.0565
Residual 1.0298 1.0148

where the estimate for the std dev of the random effects is 1.0565 while
the value used for simulating the data was 4.0.

I suspect that the real problem is with the mediocre estimation methods
used in lmer. I built the model in ADMBRE
and estimated the parameters. The estimate for sigma_u was 3.8290
4 sigma_u 3.8290e+00 1.3538e+00

Here is the TPL file. Note that I sorted the R ouput data by factor
value i.e by random effect index.


DATA_SECTION
init_int n
init_int m
int nm
!! nm=n*m;
vector y
vector x
init_matrix data(1,nm,1,3)
!! y=column(data,1);
!! x=column(data,2);
PARAMETER_SECTION
init_number a
init_number b
init_bounded_number sigma(.1,10.)
init_bounded_number sigma_u(.1,10.)
random_effects_vector u(1,m)
objective_function_value f
PROCEDURE_SECTION
int ii=0;
for (int i=1;i<=m;i++)
{
f1(a,b,u(i),sigma,sigma_u,ii);
}

SEPARABLE_FUNCTION void f1(const prevariable& a,const prevariable&
b,const prevariable& ui,const prevariable& sigma,const prevariable&
sigma_u,int& ii)

dvariable v=square(sigma);
f+=0.5*square(ui);
for (int j=1;j<=n;j++)
{
ii++;
dvariable pred=a+b*x(ii)+sigma_u*ui;
dvariable diff=y(ii)-pred;
f+=log(sigma)+0.5*square(diff)/v;
}

I wonder if the R experts are interested?






_______________________________________________
Users mailing list
Us...@admb-project.org
http://lists.admb-project.org/mailman/listinfo/users

--
You received this message because you are subscribed to the Google Groups "admb-users" group.
To post to this group, send email to admb-...@googlegroups.com.
To unsubscribe from this group, send email to admb-users+...@googlegroups.com.
For more options, visit this group at http://groups.google.com/group/admb-users?hl=en.

[H. Skaug]

Dave,


To check out your claim I compared to

lme(y ~ x, random=~1|fac,method="ML")

(the old mixed model function in R, tailor made for linear mixed models),
and got 3.829014 for the standard deviation of the random effect,
which is the same
as what you got. Note however, that these are ML estimates,
while what lmer and lme calculate by default are REML estimates.
I modified your ADMB program to calculate REML estimates and got
sigma_u=4.42137953753, which agrees with what I get from

lme(y ~ x, random=~1|fac,method="REML")

Hence, it seems that the optimization method in lmer() is not accurate enough
in this case, which is also indicated by the warning
"In mer_finalize(ans) : false convergence (8)".


I think this is a good example of why it is important to have a
good optimizer built into the package, and I will include it
in the ADMB example collection.


Hans