[ADMB Users] glmm.admb-R. Desired function, a call for calculation of 'lambda' to be used in model validation as in (fitted(Model)) vs. observed(data).

[Schweizer, Peter E.]

Dear colleagues,

We find the glmm.admb a very interesting tool for ecological modeling. Unfortunately, the at this time still rather sparse documentation of its R-application make for somehow progress in our analysis.


A few days ago we posted a question regarding the calculation and interpretation of 'residuals(Model)' in glmm.admb-R, and Hans Julius Skaug (many thanks Hans) responded. Please see below our initial question;

> > We are using glmmADMB in R to model land cover and water

> > quality influence on species diversity of fishes within a study area with

> > several subregions.

> >

> > We defined subregion as a random factor and also ask for individual

> > intercepts for the different subregions.

> >

> > A 'global' model for overdispersed count data was formulated as

> >

> > GM<-glmm.admb( N_Species ~ b1 + b2 + b3 + ...+ bn +

> Subregion, random = ~ 1,

> > group="Subregion", data=input, family="nbinom")

> >

> > We subsequently evaluated several candidate models that

> represent various

> > subsets of variables from the global model.

> >

> > Our input file is A1, with A1$NO representing the observed number of

> > species. During the process of examining model performance we used

> >

> > Observed - Predicted   (A1$NO -(fitted(best)))  for the

> 'best' model based

> > on lowest AICc to derive residuals for predicted Nspecies.

> However, using

> > 'residuals(best)' produced considerable different (smaller)

> values which we

> > find somehow puzzling. Are we wrong to assume that

> (Nspecies predicted by

> > 'best' model,  + residuals(best)) should add up to Nspecies observed

> > (A1$NO)?

> ----------------------------------------------------------------------------------------------------------
Hans Julius Skaug kindly provided the following answer;
> I think residuals(best) returns

>

> [A1$NO -(fitted(best)] / SD

>

> where SD is the standard deviation, which depends on the distribution

> at hand. The code inside glmm.admb that determines SD is:

>

>    tmpsd <- switch(family, poisson = sqrt(lambda), nbinom =

> sqrt(lambda *

>         (1 + lambda/out$alpha)), binomial = sqrt(out$fitted *

>         (1 - out$fitted)))
> ---------------------------------------------------------------------------------------
Now, since model validation for any application is an essential component of the modeling process, we are asking the ADMB community: would be possible to modify the glmm.admb R-package in the near future so that lambda can be provided in the output?



Ideally, a desired function to be developed would be a call that provides fitted(model) vs. 'Observed(model)' [=measured data from data input file, something akin to Predicted vs.Observed].



Also, at current glmm.admb-R output provides alpha as a measure of dispersion of the negative binomial distribution but without a stated lambda value, derivation of SD to calculate a P/O fit is still a challenge to be solved.



I'm sure that other colleagues in ecological research would appreciate such contribution too ...



Comments and suggestions are welcome, and thank you for your time.



Cheers,



Peter

-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.admb-project.org/pipermail/users/attachments/20110113/4b553e33/attachment.html>