[ADMB Users] ICC calculation in glmmADMB R package

Tue Sep 13 06:04:06 PDT 2011

On 09/13/2011 07:09 AM, Sen Li wrote:
> Dear Rafael
> 
> Thanks for the quick reply, it does provide me with something. But for
> the estimation of ICC, I don't know if the information is enough.
> 
> For example, when using VarCorr() in the "nlme" package, variances of
> both "intercept" and "residual" can be provided. The ICC is therefore
> calculated by: between group variance or intercept/(between group
> variance or intercept + within group variance or residual) (pp52-53,
> "Multilevel Modeling in R (ver 2.3)" by P Bliese).
> 
> When using "test2$S", I got this:
>> test2$S
>             (Intercept)
> (Intercept)      1.6267
> 
> If this gives me the "between group variance or intercept", it seems
> that I still need to get the residual variance for ICC calculation...
> How can I get this value? using var(resid(test2))?
> 
> I'm still new to R and multilevel modelling... apologise for being stupid.
> 
> With many thanks,
> Sen
> 
> 

  Not stupid at all.

  One minor comment and one major one:

  minor:  VarCorr() is available in the new (alpha) version of glmmADMB
(currently at version 0.6.4), which I hope to make the default version
soon.  But it does indeed return the same thing as test2$S in your
example above.

  major: it's not so obvious how one would compute ICC in the mixed
model context you are using -- see for example
http://r.789695.n4.nabble.com/icc-from-GLMM-td827072.html and follow
links there: Browne et al 2005 (see below) looks particularly handy.
The methods on pp. 603-605 there are derived for a binomial rather than
a negative binomial model; one could probably derive the equivalents for
a NB model (or, more closely, a lognormal-Poisson model which you would
get by fitting a Poisson with an additional observation-level random
effect) -- it would take me a few hours (give or take), might take you
longer ... Nakagawa and Schielzeth 2010 discuss measures of
repeatability and heritability in non-Gaussian models, which are related
to the ICC.

  Based on a little more poking around (the Google results from "ICC
'negative binomial' are reasonably enlightening) I think the simplest
(??) way to recover an ICC for overdispersed count data would be to fit
a lognormal-Poisson model and then take the ICC as

    group variance/(group variance + individual-level variance + 1)

where the 1 would represent the Poisson variance.

  Bottom line: unfortunately, I haven't found a simple, universal answer
to this.

@article{browne_variance_2005,
	title = {Variance partitioning in multilevel logistic models that
exhibit overdispersion},
	volume = {168},
	issn = {{1467-985X}},
	url =
{http://onlinelibrary.wiley.com/doi/10.1111/j.1467-985X.2004.00365.x/abstract},
	doi = {10.1111/j.1467-985X.2004.00365.x},
	abstract = {Summary.  A common application of multilevel models is to
apportion the variance in the response according to the different levels
of the data. Whereas partitioning variances is straightforward in models
with a continuous response variable with a normal error distribution at
each level, the extension of this partitioning to models with binary
responses or to proportions or counts is less obvious. We describe
methodology due to Goldstein and co-workers for apportioning variance
that is attributable to higher levels in multilevel binomial logistic
models. This partitioning they referred to as the variance partition
coefficient. We consider extending the variance partition coefficient
concept to data sets when the response is a proportion and where the
binomial assumption may not be appropriate owing to overdispersion in
the response variable. Using the literacy data from the 1991 Indian
census we estimate simple and complex variance partition coefficients at
multiple levels of geography in models with significant overdispersion
and thereby establish the relative importance of different geographic
levels that influence educational disparities in India.},
	number = {3},
	journal = {Journal of the Royal Statistical Society: Series A
{(Statistics} in Society)},
	author = {Browne, W. J and Subramanian, S. V and Jones, K. and
Goldstein, H.},
	month = jul,
	year = {2005},
	keywords = {Contextual variation, Illiteracy, India, Multilevel
modelling, Multiple spatial levels, overdispersion, Variance partition
coefficient},
	pages = {599--613}
}

@article{nakagawa_repeatability_2010,
	title = {Repeatability for Gaussian and {non-Gaussian} data: a
practical guide for biologists},
	issn = {14647931, {1469185X}},
	shorttitle = {Repeatability for Gaussian and {non-Gaussian} data},
	url = {http://doi.wiley.com/10.1111/j.1469-185X.2010.00141.x},
	doi = {10.1111/j.1469-185X.2010.00141.x},
	journal = {Biological Reviews},
	author = {Nakagawa, Shinichi and Schielzeth, Holger},
	month = jun,
	year = {2010},
	pages = {no--no}
}

> On 13/09/2011 12:20, Rafael Mares wrote:
>> Is test2$S what you are looking for? I think it is similar to (or the
>> same as?) VarCorr() but for glmm.ADMB models.
>>
>>
>> -- 
>> Rafael Mares
>> Large Animal Research Group (LARG)
>> Department of Zoology
>> University of Cambridge
>> Downing Street
>> Cambridge
>> CB2 3EJ
>>
>>
>>
>> On 13 September 2011 09:33, Sen Li<sen.li at uclouvain.be>  wrote:
>>> Dear Ben,
>>>
>>> Thanks for the reply. I tried VarCorr() before, but it turned out to be
>>> inapplicable with gimm.ADMB objects.
>>>
>>> Here follows the outputs of my model, of VarCorr(), and of
>>> sessionInfo():
>>>
>>>> print(test2,sd_S_print=T)
>>> GLMM's in R powered by AD Model Builder:
>>>
>>>    Family: nbinom
>>>    alpha = 1.0808
>>>
>>> Fixed effects:
>>>    Log-likelihood: -690.741
>>>    AIC: 1395.482
>>>    Formula: y ~ x1 + x2 + x3
>>>   (Intercept)           x1                x2                    x3
>>>     11.048000    -0.025435   -45.852000    -0.061126
>>>
>>> Random effects:
>>>    Grouping factor: Scode
>>>    Formula: ~1
>>> Structure: Diagonal matrix
>>> (Intercept)
>>>     1.275421
>>>
>>> Covariance matrix of random effects vector (left) and corresponding
>>> standard
>>> deviations (right):
>>>
>>>              (Intercept)     (Intercept)
>>> (Intercept)      1.6267         0.77551
>>>
>>> Note: The diagonal elements of the above left matrix are variances, NOT
>>> std's.
>>>
>>> Number of Observations: 125
>>> Number of Groups: 51
>>>> VarCorr(test2)
>>> Erreur dans UseMethod("VarCorr") :
>>>    pas de méthode pour 'VarCorr' applicable pour un objet de classe
>>> "glmm.admb" ## the method is not applicable to an object of class
>>> "glmm.admb"##
>>>> sessionInfo()
>>> R version 2.13.0 (2011-04-13)
>>> Platform: x86_64-pc-mingw32/x64 (64-bit)
>>>
>>> locale:
>>> [1] LC_COLLATE=French_Belgium.1252  LC_CTYPE=French_Belgium.1252
>>> [3] LC_MONETARY=French_Belgium.1252 LC_NUMERIC=C
>>> [5] LC_TIME=French_Belgium.1252
>>>
>>> attached base packages:
>>> [1] stats     graphics  grDevices utils     datasets  methods   base
>>>
>>> other attached packages:
>>> [1] nlme_3.1-101   glmmADMB_0.5-2 MASS_7.3-12
>>>
>>> loaded via a namespace (and not attached):
>>> [1] grid_2.13.0     lattice_0.19-23 tools_2.13.0
>>> I'd appreciate any further help.
>>>
>>> Best,
>>> Sen
>>>
>>> _______________________________________________
>>> Users mailing list
>>> Users at admb-project.org
>>> http://lists.admb-project.org/mailman/listinfo/users
>>>
>>>
> 
> 
> _______________________________________________
> Users mailing list
> Users at admb-project.org
> http://lists.admb-project.org/mailman/listinfo/users