problem met in vol1 example Leuk
Jiqiang Guo
Exception: Error in function stan::prob::normal_log<d>(d): Random
variate x is -nan, but must be finite!
Diagnostic information:
Throw location unknown (consider using BOOST_THROW_EXCEPTION)
Dynamic exception type:
boost::exception_detail::clone_impl<boost::exception_detail::error_info_injector<std::domain_error>
>
std::exception::what: Error in function stan::prob::normal_log<d>(d):
Random variate x is -nan, but must be finite!
I could not figure out what is wrong with the model specification. So
help needed.
--
Jiqiang
Daniel Lee
try to help hunt it down.
Daniel Lee
nan. I'm trying to figure out how it gets set to nan.
Bob Carpenter
>> Error message during running:
>>
>> Exception: Error in function stan::prob::normal_log<d>(d): Random
>> variate x is -nan, but must be finite!
This is exciting -- we're printing reasonable
messages for errors already.
Given that the only normal in the model is
beta ~ normal(0, 10);
and the beta isn't declared to be transformed, the beta
parameter had to come in as a NaN value. I don't see how
this could come about -- the only thing that happens to it
is adding momentum amounts and the gradient of normal_log(0,10)
shouldn't be too sharp.
Is it repeatable, Daniel? It always goes to NaN
in the second step?
One issue is whether we then want to just catch
these exceptions in the sampler loop, report
them, then try again (up to some limit so we
don't go into an infinite loop).
I'm not sure what JAGS does, but I've often seen
pages and pages of error messages about NaN values
found in evaluations (usually in gamma functions).
- Bob
Daniel Lee
Yes, it is repeatable.
I think we may want to catch it at some point in the future, but for
the time being, it's useful to fail fast, I think.
Daniel
Daniel Lee
We haven't bounds checked the Poisson distribution yet. The rate
parameter, lambda, is 0 at some point.
Jiqiang Guo
Bob Carpenter
all these boundary conditions that might pop up.
We can set
poisson_log(N,lambda)
= 0.0 if N = 0, lambda = 0.0
= error if N > 0, lambda = 0.0
= error if lambda < 0 or N < 0
= business as usual otherwise
- Bob
Jiqiang Guo
Daniel Lee
been updated.
The model seems to be running now. Jiqiang, will you let me know if it
does the right thing?
Daniel
Jiqiang Guo
Bob Carpenter
> > dpois(1, lambda = 0)
> [1] 0
> > dpois(0, lambda = 0)
> [1] 1
> > dpois(2, lambda = 0)
> [1] 0
> > dpois(-2, lambda = 0)
> [1] 0
That makes sense mathematically if you
think of Poisson(N|lambda) being defined over
all integers for all N no matter what lambda is.
I think Boost has the same kind of behavior.
Also, we're following Mssrs. Gelman, Carlin, Stern,
and Rubin, the appendix of whose book states
that Poisson(N|lambda) is defined for lambda > 0 (strictly)
and N = 0, 1, 2, ...
The problem is that even if lambda is a parameter
which we bound > 0, we can get 0 due to underflow
when we do the exp(x) reverse transform for
x < -400 or so.
We also need gradients. The problem is
that poisson_log(N,lambda) doesn't have a gradient
if N is negative, or if lambda is not greater than
zero.
We also want to be more helpful in flagging client
errors. If they put in count data with negative
counts, we want to trap it as early as possible as
an error.
Another issue for us is that we want to be able
to
Like Boost, R handles the RHS differently.
> dpois(-2,-2)
[1] NaN
Warning message:
In dpois(x, lambda, log) : NaNs produced
Unlike Boost, R would make Matt happy, because it
treats the location and sampled value symmetrically.
Boost complains if you do y ~ norm(mu,sigma) with
anything but a finite mu, though it'll allow infinite y.
(That makes sense given their design because they
construct an object for norm(mu,sigma) and catch the
error at construction time if you've used illegal
parameters.)
> So it might be better that when lambda=0,
> poisson_log(N, lambda) gives log(0) when N is not 0.
> I am not sure about this since even I am not sure about
> log(0) is legitimate.
log(0) produces -infinity in IEEE arithmetic.
We'd like to make sure our distribution code handles
the infinities well, too. For instance, R also does
sensible things with all the infinity cases for the normal:
> dnorm(1/0,0,1)
[1] 0
>
> dnorm(1/0,0,1,log=T)
[1] -Inf
> dnorm(0,1/0,1)
[1] 0
> dnorm(0,1/0,1,log=T)
[1] -Inf
> dnorm(0,0,1/0)
[1] 0
> dnorm(0,0,0/1)
[1] Inf
- Bob
Daniel Lee
lambda == 0 (regardless of N) is a domain error?
or is it
lamdba == 0 when N > 0 is a domain error?
Daniel
Bob Carpenter
> I think at the end of that, we've decided that:
> lambda == 0 (regardless of N) is a domain error?
> or is it
> lamdba == 0 when N> 0 is a domain error?
Here's my thinking:
log_poisson(n,lambda)
= log Poisson(n|lambda) if +infinity > lambda > 0 and n >= 0
= domain error if lambda < 0 or n < 0
= domain error if n is NaN or lambda is NaN
The NaN error and other errors may just lead to
NaN propagation with the right traits set.
That leaves a few boundary conditions undecided.
= domain error if lambda = +infinity
If it were zero, the partial derivative w.r.t. lambda is
-infinity, and when we do +infinity + -infinity, we'll
get a NaN. So I think it's safe to leave this
an error.
= 0 if n = 0, lambda = 0
Just branching and returning 0 makes the gradient 0,
too. It could also throw an error.
= domain error if n > 0, lambda = 0
This could return -infinity if the gradient
makes sense:
lim(L -> 0) d/d.lambda log Poisson(n|lambda)
I didn't work through whether this limit makes
sense or blows up. The gradient would probably
have to be plumbed in by hand.
- Bob
Wei Wang
Wei
Matt Hoffman
few scenarios are
lambda=0, n=0: logp=0, grad_lambda=–1 (because 0log0 is basically 0)
lambda=0, n>0: logp=-inf, grad_lambda=inf (for one thing, because
log(1-1)=-inf and log(1-1+epsilon)>-inf)
lambda=inf, n<inf: logp=-inf, grad_lambda=-inf (because
log(gamma(lambda)) grows much faster than n*log(lambda))
lambda=inf, n=inf: domain error
Matt
Bob Carpenter
a whole bunch of cdf gradients.
I didn't work through the actual limits before, but
my intuition was that poisson(n|lambda) would converge
to a delta function as n->0 and lambda->0.
0 log 0 is another case we might want to build
into our auto-diff. There's no way to get this
out of basic log evaluation, though, because log(0) is
negative infinity and you evaluate compositionally.
Now if we had expression templates we might be able
to make it fly, but even then, it's a semantic not
a syntactic calc for the value of x in log(x) unless
the x is a structural 0 in a model somewhere.
- Bob
On 12/12/11 2:35 AM, Matt Hoffman wrote:
> The limit values that the Poisson log-pdf actually gives in those last
> few scenarios are
> lambda=0, n=0: logp=0, grad_lambda=�1 (because 0log0 is basically 0)