propto density summands & unit testing

[Bob Carpenter]

I'm in the middle of fixing multi normal and
noticed that there are not enough tests for propto.

You can drop a summand in a log probability if and
only if both

1) propto == true, and
2) there are no stan::agrad::var variables in the
summand

Here's the first summand in the definition of multi_normal:

if (!propto)
lp += NEG_LOG_SQRT_TWO_PI * y.rows()
- 0.5 * log(Sigma.determinant());

This needs to be broken out into:

if (!propto)
lp += NEG_LOG_SQRT_TWO_PI * y.rows();

if (!propto || !is_constant<T_covar>::value)
lp -= 0.5 * log(Sigma.determinant());

I'll add the is_constant<> metaprogramming for
vectors and matrices to support the static tests.

The other conditions need to be refined in the
same way. If we include all the summands, we'll
get the right answer -- it'll only cost us some speed
and memory. So feel free to leave things without
the propto calc -- we can optimize later in a way
that'll be backward compatible.

We will need tests for propto. I think these should
be in the form of log difference tests, exploiting
the fact that we've tested foo<false> (full calc),
and that:

foo<true>(theta1) - foo<true>(theta2)

== foo<false>(theta1) - foo<false>(theta1)

So we can choose some pairs of theta1 and theta2,
then work through all the permutations of var
and constant testing the above.

This'll avoid the need to test exact values of
foo<true> (proportional calc), which may change over
time as we optimize. (It won't if we minimize
the calcs the first time through.)

And the whole reason I dove into the class in the
first place is that we need to replace the direct calls
to Eigen's matrix ops with our own wrappers to take
care of promotion to auto-dif variables (alternatively,
we could promote within the probability calcs, but that
seems less modular).

- Bob

[Daniel Lee]

My fault on the multi_normal.

I realized I made some errors -- I'll go back through and fix them.

Regarding tests, I've started to test the distributions propto in
src/test/agrad/distributions_*_test.cpp.

[Daniel Lee]

Could you take a look at src/test/agrad/distributions_normal.cpp?

I've used a test fixture (TEST_F). Do you think something like this
will work? Should we trust the auto-dif completely or should we check
that too?

Daniel

PS. If you get to my mistakes before I get back through them, please
feel free to update them.

[Daniel Lee]

I got through taking another look at all the distributions. There is
still a lot of testing to be done to ensure that everything is
correct, but I'm a lot more confident in it now.

Something I'm not quite clear about is how we deal with the random
variable, y, in the propto calculations. If y is a double, which I
expect it to be often enough, should that not be used in the propto
calculation? Or should we just use it in our calculations all the
time?

Daniel

[Bob Carpenter]

I think it's OK to have unit tests depend on lower-level
functions that are themselves tested. We don't need to
test compound functionality the pieces of which are unit
tested. So just testing that we've defined the probability
densities so they can be auto-diffed is enough -- we don't
need to check the results (of course, it couldn't hurt and
redundant testing with different approaches is in general
a good thing).

There's nothing special about the y in y ~ Foo(theta),
You'll have to include it when it touches the params.
So it's usually still there following the same reules.

For instance, for log Beta(y|alpha,beta):

summand include if variable
1) - log B(alpha,beta) alpha | beta
2) + (alpha - 1) log(y) alpha | y
3) + (beta - 1) log(1 - y) beta | y

Here (1) is a pure normalizing term only involving the
distribution arguments, but y touches the parameters
in terms (2) and (3). The big savings is if alpha and
beta are constant, because (1) involves log gamma functions
and derivatives digamma functions as well.

- Bob

PS: I'm on my way into the office.

[Daniel Lee]

Got it. For the most part, I have been doing it that way. I'll just
make sure everything is done that way.

I realized my tests don't quite work the way I want. I've only been
testing proportionality by comparing the log diff for two y's, but
what I need to do is to look at proportionality for the parameters
that are vars.

Daniel

[Bob Carpenter]

On 12/19/11 12:51 PM, Daniel Lee wrote:
> Got it. For the most part, I have been doing it that way. I'll just
> make sure everything is done that way.
>
> I realized my tests don't quite work the way I want. I've only been
> testing proportionality by comparing the log diff for two y's, but
> what I need to do is to look at proportionality for the parameters
> that are vars.

If we're being careful, we'll look at all the permutations.
I'll go in and do that for the normal as an example. It's
going to be a lot of work.

We only really need to do it for ones where we've optimized
so the propto works differently than the non-propto version.

- Bob

[Daniel Lee]

Part of me wants to have all the propto tests in just so it's there as
a failsafe as we start optimizing. Think it's worth it?

[Bob Carpenter]

It's definitely worth it -- just a matter of how
much stamina you have to push it through.

I've organized the normal tests so there's not so
much code involved. I'll ping you when it's checked
in.

- Bob