Dear List,
I have a related question/comment about DIC but it is not specific to
Winbugs. I am using the mcmc component of the AD model builder software to
estimate posterior distributions for a range of hierarchical models and then
computing DIC statistics for each model myself (in R). When computing Dhat I
use the most likely (lowest deviance) estimate from the posterior sample.
This seems to work, as the ranking of converged models via DIC fits my
expecations, pD is always less than the total number of parametes for each
model (K), and the difference between pD and K for models with equal K lines
up with differences in the structural rigidity of the models.
From the documentation in WinBugs, Dhat is computed from the mean parameter
values of the posterior, slightly different than my approach of using the
lowest deviance from the posterior. Gelman et al. (2nd edition, p. 181) are
a bit vague on how Dhat should computed:
"One could define Dhat using a point estimate for theta such as the mean of
the posterior simulations" (sic)
By using the word 'could', are they implying that other statistics can be
used? Is there some theoretical reason why the mean parameter values of the
posterior must be used to compute Dhat or should Dhat just reflect a point
estimate of the lowest deviance?
Any help/comments would be much appreciated.
Cheers,
Josh Korman, MSc
PhD candidate, Department of Zoology, University of British Columbia
3560 W 22nd Ave.
Vancouver, BC
V6S 1J3
Canada
tel. (604) 737-8314
fax. (604) 737-8324
email: [log in to unmask]
----- Original Message -----
From: "Jonathan Rhodes" <[log in to unmask]>
To: <[log in to unmask]>
Sent: Wednesday, January 17, 2007 6:14 PM
Subject: [BUGS] DIC and DBar (Summary of Responses)
Dear List,
I recently posted the following question regarding Dbar in the calculation
of DIC, to which I received several responses.
"I am using DIC to compare alternative (nested) generalised linear
mixed-effects models, with the alternative models containing different
fixed-effects (the random-effects do not vary between models). However, I
notice that, for example, when comparing two models, with one model having
one less fixed-effects parameter than the other (the model with fewer
parameters being nested within the model with more parameters) that often
Dbar is higher for the model with more parameters than for the model with
fewer parameters. Intuitively I would expect Dbar to be lower for the model
with more parameters. Am I missing something here, or is this a sign of poor
convergence or mixing (some of my parameters have reasonably high
autocorrelations)? Has anyone else come across this problem?"
From the responses I received this seems to be a fairly common problem.
Quite a number of respondents suggested that it is likely to be a
convergence problem, as I originally thought it might be. Having had a
closer look at my models (they are mostly logistic regression with normal
random-effects) this certainly seemed to be true in some cases (e.g. due to
inappropriate random-effects, or high collinearity between covariates). For
other models convergence appeared to be quite good, but I still sometimes
got an increase in Dbar when adding a parameter. David Spiegelhalter
suggested that, for random-effects models, it might theoretically be
possible for Dbar to increase. Also, Thomas Jagger provided a simple
example, using a normal distribution, showing that Dbar can theoretically
increase when a parameter is added, but that Dhat must always decline (at
least for the model used in his example). His argument appears to make sense
and for my models that converged well, this is exactly what I see: sometimes
I get a slight increase in Dbar, but Dhat always declines. So, it does not
seem entirely unreasonable for Dbar to increase when a parameter is added to
the model.
A number of other respondents (Murray Aitken and Seth Wenger) suggested
alternative model selection approaches using marginal likelihoods or
cross-validation that people may be interested in.
I have pasted a selection of the most relevant responses below.
SELECTED RESPONSES
From: Vallejo, Roger [[log in to unmask]]
Sent: Thursday, 14 December 2006 2:13 AM
To: Rhodes, Jonathan (CMAR, Hobart)
Subject: RE: [BUGS] DIC and Dbar
This is typical sign of poor mixing or convergence. The reasons could be
several. Autocorrelations as you indicate is one, rate of thinning going
from 1/10 to 1/50 can help and increase the number of runs, poor initial
parameter values, etc.
Roger
-----
From: David Spiegelhalter [[log in to unmask]]
Sent: Monday, 18 December 2006 9:49 PM
To: Rhodes, Jonathan (CMAR, Hobart)
Subject: Re: DIC and Dbar
Jonathan
With random effects models, I suppose adding fixed effects could
theoretically lead to a higher Dbar, as the random effects estimates will
also be changed and somehow fit the data worse with more fixed effects.
However I think such behaviour would be very unusual, and maybe suggests
something odd with the model/data?
David
-----
From: Thomas Jagger [[log in to unmask]]
Sent: Wednesday, 20 December 2006 6:24 PM
To: Rhodes, Jonathan (CMAR, Hobart)
Subject: RE: [BUGS] DIC and Dbar
Hi Jonathan,
You may already have the answer, but consider that the DIC= 2*Dbar - Dhat
This implies that if the increase in Dbar is less than 1/2 the decrease in
Dhat, that in fact, the new model fits better, even though Dbar has
increased. Dhat should always decrease, but Dbar does not need to.
Consider even the simple example, suppose we have N samples (x) of a normal
distribution with unknown mean and variance=1. Assume we have no parameters,
and we assume that the mean is zero, then Dhat=Dbar, and deviance (removing
the log(2*pi) is sum(x^2). Now suppose we assume a mean parameter, mu, with
non-informative prior then Dhat = sum((x-mean(x))^2)= sum(x)^2- N*mean(x)^2.
The distribution of the posterior mean has a normal distribution with
mean=mean(x) and variance=1/N.
Each term in the likelihood has the form (x-mu)^2, the posterior mean of
(x-mu)^2 is then
x^2-E(mu)*x+E(mu^2)=x^2-mean(x)*x-mean(x)^2+1/N.
When you sum all the terms up you get:
Dbar=sum(x^2)-2*N*mean(x)^2+N*mean(x)^2+1/N=sum(x^2)-N*mean(x)^2+1=Dhat+1
which makes sense as Dbar-Dhat=Pd=1.
If the mean(x)=0, which might be reasonable, then dbar increases by 1.
The expectations over all x for N iid N(0,1) samples, give Dhat=N-1 and
Dbar=N and DIC=N-1+2=1.
Now the change in Dbar is 1-N*mean(x)^2 = 1-Y where under the assumption
that x has N iid N(0,1) random variables, then Y has a chi-square
distribution with one degree of freedom so Dbar has mean of 0 and variance
of 2. Note that the change in DIC=2-Y, so that the DIC can increase by at
most 2 (Dhat stays the same).
Hope this helps.
Tom Jagger
-----
From: Seth Wenger [[log in to unmask]]
Sent: Thursday, 14 December 2006 4:06 AM
To: Rhodes, Jonathan (CMAR, Hobart)
Subject: Re: DIC and Dbar
Jonathan,
I have experienced this as well, and have talked to others who have reported
this. I'm no statistician, but I'm going to give you my current thinking on
it, which is quite possibly incorrect.
- It seems that WinBUGS calculates the likelihood with the random effects
included, effectively as if they were fixed effects. This differs from
traditional approaches, in which the random effects are integrated out of
the likelihood.
- DIC is an attempt to address this, but doesn't quite solve it, as you have
observed.
- Possible solutions are (a) to create your own likelihood node, but I'm
unclear on the most appropriate way to do this, or (b) to use
cross-validation predictive success as the model selection criterion; the
REs are excluded from the predictions. I have done the latter and I think
it is entirely appropriate. But it is a fair bit of work, so for practical
purposes you'll need to limit it to 3-fold CV or something similar;
leave-one-out will probably not be an option.
- The entire issue makes me cautious about incorporating random effects, and
I try to think very carefully about why I'm including random effects in the
model. Is it a tool for minimizing bias associated with spatial
correlations, for example? I think it's worth giving serious thought.
I'd appreciate you sharing any good responses you receive from folks who
understand this issue better than I do.
Thanks!
Seth Wenger
-----
From: Murray Aitkin [[log in to unmask]]
Sent: Thursday, 14 December 2006 11:50 AM
To: Rhodes, Jonathan (CMAR, Hobart)
Cc: [log in to unmask]; Charles Liu; Tom Chadwick
Subject: Re: [BUGS] DIC and Dbar
Jonathan - There's an alternative to DIC which you might investigate - it
uses directly the posterior distribution of the deviance (not the complete
data deviance) without penalty, which isn't needed. The theory is in Aitkin,
Boys and Chadwick Statistics and Computing (2005) 15, 217-230. I'll send you
separately a paper by Aitkin, Liu and Chadwick on comparing models for
two-level data using this approach.
Cheers,
Murray
-----
Cheers,
Jonathan Rhodes
Postdoctoral Fellow
CSIRO Marine and Atmospheric Research
Castray Esplanade
Hobart
TAS 7000
Australia
Tel: +61-(0)3-62325113
Fax: +61-(0)3-62325000
Postal address:
CSIRO Marine and Atmospheric Research
GPO Box 1538
Hobart
TAS 7001
Australia
Links: www.cmar.csiro.au (CSIRO Marine and Atmospheric Research)
-------------------------------------------------------------------
This list is for discussion of modelling issues and the BUGS software.
For help with crashes and error messages, first mail [log in to unmask]
To mail the BUGS list, mail to [log in to unmask]
Before mailing, please check the archive at
www.jiscmail.ac.uk/lists/bugs.html
Please do not mail attachments to the list.
To leave the BUGS list, send LEAVE BUGS to [log in to unmask]
If this fails, mail [log in to unmask], NOT the whole list
-------------------------------------------------------------------
This list is for discussion of modelling issues and the BUGS software.
For help with crashes and error messages, first mail [log in to unmask]
To mail the BUGS list, mail to [log in to unmask]
Before mailing, please check the archive at www.jiscmail.ac.uk/lists/bugs.html
Please do not mail attachments to the list.
To leave the BUGS list, send LEAVE BUGS to [log in to unmask]
If this fails, mail [log in to unmask], NOT the whole list
|
