summing out int parameters
Bob Carpenter
parameters wherever we can. It'll be faster to execute per
sample and convergence/exploration should be faster,
too.
(Is this what people call Rao-Blackwellization? I've also
seen it called "collapsed" sampling in the LDA setting
when you integrate out continuous params.)
An example is the Hearts model (BUGS Vol2 Example)
that Jiqiang just checked in (I just love the
diagram for GIT checkins on
http://code.google.com/p/stan/source/list
# Hearts
data {
int(0,) N;
int(0,) x[N];
int(0,) y[N];
int(0,) t[N];
}
parameters {
real alpha;
real delta;
int(0,) state[N]; // these are discrete parameters. the support is 0 or 1 actually. ??
}
transformed parameters {
real(0,1) theta;
theta <- inv_logit(delta);
}
model {
real yap[2];
yap[1] <- inv_logit(alpha);
yap[2] <- 0;
alpha ~ normal(0, 100);
delta ~ normal(0, 100);
for (i in 1:N) {
state[i] ~ bernoulli(theta); // either 0 or 1
y[i] ~ binomial(t[i], yap[state[i] + 1]);
}
}
(Note: the int(0,) should be int(0,1) for state[N] because
it's binary; but that's the parameter we're summing out
so we don't need it at all.)
What we need to do is replace the inner loop statements:
state[i] ~ bernoulli(theta); // either 0 or 1
y[i] ~ binomial(t[i], yap[state[i] + 1]);
with the version that sums out state[i]:
lp__ <- lp__
+ log_sum_exp(log(theta) + binomial_log(y[i],t[i],yap[1]),
log(1.0 - theta) + binomial_log(y[i],t[i],yap[2]));
The increment expression is just an arithmetically stable
version of:
log(theta * binomial(y[i],t[i],yap[1])
+ (1-theta) * binomial(y[i],t[i],yap[2]);
I haven't written this section of the manual yet. But
it's also on the to-do list. You can leave this example
to me and I'll make sure it works when I write the manual
section.
- Bob
Jiqiang Guo
Andrew Gelman
A
Bob Carpenter
On 4/12/12 10:28 PM, Andrew Gelman wrote:
> I think that solving all the models in the Jags manual is a good start,
There are no models in the JAGS manual per se.
They released a code tarball that has translations
of the "BUGS models" (see previous message) to JAGS,
but there's no documentation, just the code.
> also solving all the models we care about.
Is there a list somewhere? Are there models we care
about that can't be expressed in BUGS/JAGS?
We should be able to handle anything BUGS/JAGS
does and quite a bit more in terms of expressiveness.
> Probably best to solve models one at a time, developing general tools as needed, rather than trying to anticipate
uses ahead of time.
Too late for that given that our initial menu was
already all of the "BUGS models" and "everything we
care about" (e.g., speed dating, folates, etc.).
Also, it's not a very good strategy
for developing a programming language -- you wind up
with a Frankenstein monster like Perl if you don't
plan ahead for how things are going to be extended.
The big issue on the expressiveness front is Bayesian
non-parametric models. We haven't thought much at
all about those either for expression or implementation.
- Bob
Jiqiang Guo
> also solving all the models we care about.
Is there a list somewhere? Are there models we care
about that can't be expressed in BUGS/JAGS?
We should be able to handle anything BUGS/JAGS
does and quite a bit more in terms of expressiveness.
Bob Carpenter
On 4/13/12 1:53 PM, Jiqiang Guo wrote:
>
>
> On Fri, Apr 13, 2012 at 1:37 PM, Bob Carpenter <ca...@alias-i.com <mailto:ca...@alias-i.com>> wrote:
> Openbugs's official example webpage:
> http://www.openbugs.info/w/Examples
That's the same HTML as on the WinBUGS site.
It's what I think of as the "official" list
of BUGS models (in both senses of model).
> JAGS does not have an official stuff for examples, I suppose.
> But http://sourceforge.net/projects/mcmc-jags/files/Examples/2.x/
> has some BUGS vol1 and vol2 examples.
Right -- that's what I was referring to.
> And JAGS' author has a R package for BUGS examples:
> https://r-forge.r-project.org/projects/bugs-examples/
> This package includes more examples than the previous one, but
> it does not explicitly say which examples are not supported in JAGS.
Thanks -- I hadn't seen that.
> This package is written in a way you could run any example using
> either JAGS or BUGS. And I guess that is how the speed
> comparison is done (post this link again:
> http://martynplummer.wordpress.com/2010/09/20/how-fast-is-jags/).
I should look a bit closer at OpenBUGS. I just went
straight to using JAGS from WinBUGS. If they're testing
from R, that could be a big factor in timing itself.
I'm going to break out a timing thread separately, as it's
something we need to do in the near term.
- Bob
Ben Goodrich
And JAGS' author has a R package for BUGS examples:This package includes more examples than the previous one, butit does not explicitly say which examples are not supported in JAGS.
He cached the supported models in a data.frame but I don't know how up-to-date that is. -- Ben
> require(BUGSExamples); print(bugsex) name WinBUGS OpenBUGS JAGS 1 Air TRUE TRUE TRUE 2 Alligators TRUE TRUE TRUE 3 Asia TRUE TRUE TRUE 4 Beetles TRUE TRUE TRUE 5 Biopsies TRUE TRUE TRUE 6 BiRats TRUE TRUE TRUE 7 Blockers TRUE TRUE TRUE 8 Bones TRUE TRUE TRUE 9 Camel TRUE TRUE FALSE 10 Cervix TRUE TRUE TRUE 11 Dogs TRUE TRUE TRUE 12 Dugongs TRUE TRUE TRUE 13 Dyes TRUE TRUE TRUE 14 Elevation TRUE TRUE FALSE 15 Endo TRUE TRUE TRUE 16 Epil TRUE TRUE TRUE 17 Equiv TRUE TRUE TRUE 18 Eyes TRUE TRUE TRUE 19 Eyetracking TRUE TRUE TRUE 20 Forest TRUE TRUE FALSE 21 Hearts TRUE TRUE TRUE 22 HepatitisME TRUE TRUE TRUE 23 Hepatitis TRUE TRUE TRUE 24 Hips2 TRUE TRUE TRUE 25 Hips3 TRUE TRUE TRUE 26 Hips4 FALSE FALSE FALSE 27 Huddersfield FALSE FALSE FALSE 28 Ice TRUE TRUE FALSE 29 Inhalers FALSE FALSE FALSE 30 Jama TRUE TRUE FALSE 31 Jaws TRUE TRUE TRUE 32 Kidney TRUE TRUE FALSE 33 Leukfr TRUE TRUE TRUE 34 Leuk TRUE TRUE TRUE 35 LHA TRUE TRUE FALSE 36 Line TRUE TRUE TRUE 37 Lsat TRUE TRUE TRUE 38 Magnesium TRUE TRUE TRUE 39 Mice TRUE TRUE TRUE 40 MVCAR TRUE TRUE FALSE 41 Otrees TRUE TRUE TRUE 42 OtreesMVN TRUE TRUE TRUE 43 Oxford TRUE TRUE TRUE 44 Pigs TRUE TRUE TRUE 45 Pigweights TRUE TRUE TRUE 46 Pines FALSE FALSE FALSE 47 Pollution TRUE TRUE FALSE 48 Pumps TRUE TRUE TRUE 49 Rats TRUE TRUE TRUE 50 Salm TRUE TRUE TRUE 51 Schools TRUE TRUE TRUE 52 Scotland1 TRUE TRUE FALSE 53 Scotland TRUE TRUE FALSE 54 Seeds TRUE TRUE TRUE 55 Shared TRUE TRUE FALSE 56 Sixcomp FALSE FALSE FALSE 57 Stacks TRUE TRUE TRUE 58 Stagnantc TRUE TRUE TRUE 59 StVeit TRUE TRUE FALSE 60 Surgical TRUE TRUE TRUE 61 Surgicalrand TRUE TRUE TRUE 62 t-df FALSE FALSE FALSE