added MCMC posterior_predictive_branch sampler

[danielturek]

@perrydv @paciorek Added checking for jointly posterior predictive networks of nodes and a new posterior_predictive_branch sampler to be assigned to them.

I've tested this (what I believe to be) pretty well, although admittedly no new testing was added. But I believe it's sound. Code review or testing is welcomed.

[danielturek]

@perrydv @paciorek I've reviewed all the failures carefully, because I was curious if the sampler assignment logic held up.

The first failures are cases where the new posterior_predictive_branch should be assigned, but in the tests, we're checking (usually) for the presence of either RW or conjugate samplers, to make sure the model inspection and sampler assignment works.

In test-trunc.R, in the test labelled "Test that MCMC with truncation avoids conjugate samplers". I encourage anyone interested to look at this test, and understand what's going on, with the new assignment of the posterior_predictive_branch sampler. Best fix might be making all the y terms into data.

Next failure is in test-dynamicIndexing.R. There are two tests, labelled "Testing multivariate normal-Wishart dependency conjugacy detection with dynamic indexing", and all the failures are taking place there. Maybe, a better fix for this one, would be just specifying the y as being "data", in both of the tests? Then, we'd still check for (and should assign) conjugate samplers, which is what the test is trying to test for. As is, with no data in the models, these nodes get the new posterior_predictive_branch sampler.

Finally, there are a bunch of failures again the "gold file" of exact MCMC results that I'm seeing, I think all from test-mcmc.R. I haven't looked carefully at all of these, but of course it makes sense that we'd be getting different samples now.

Thought about how to proceed?

[danielturek]

@perrydv @paciorek I fixed the testing failures.

• In tests that check for certain sampler assignments, I added "data" to the models, so the old checking for particular conjugate or non-conjugate relationships still applies (rather than the new posterior_predictive_branch sampler).
• In tests checking for exact sample results, or matching the MCMC gold file, I added a new system option MCMCusePosteriorPredictiveBranchSampler which disables use of this new sampler. This option is turned off for the test-mcmc.R tests, and it now matches old results.

I'd welcome any review on this. @perrydv also your careful eye for efficiency.

[danielturek]

@perrydv I'm going to modify this PR to use node graph IDs instead.

[danielturek]

@perrydv Unfortunately, some quick efficiency tests have shown this approach to be unusable. Back to the drawing board.

I have an alternate idea for a "bottom-up" approach to the problem, using getParents.

[danielturek]

@perrydv This checking and assignment is now "fast", and doesn't bog down the MCMC configuration process.

I'd welcome any code review, or to merge this into devel.

[paciorek]

Many of our options are framed as "actions", i.e., they have a verb in them. If "Sample" here is a verb, then it should be "jointlySample" not "jointSample". Obviously minor, but my two cents.

[danielturek]

@paciorek I added your idea of avoiding any checking, when there's no possibility whatsoever of any posterior branches. Indeed, in some cases this totally alleviates the checking time. Good idea.

[perrydv]

@danielturek Here are some suggestions on the code:

1. It would be nice to do model$getDependencies(thisCandNodeID, dataOnly = TRUE, downstream = TRUE, returnType = 'ids') once and re-use the results. It looks like it is done up to five times redundantly. One of those times has self = FALSE, which would need to be dealt with.
2. Could the three for loops be combined (which would help on the previous point)? It looks like the indToKeep and nodeIDs steps could be done in the first for loop.
3. I think the direct use of the maps internals could be replaced with model$expandNodeNamesFromGraphIDs.

[perrydv]

@danielturek Here is a possibly challenging benchmark case:

mu ~ dnorm(0, 1) # prior
y ~ dnorm(mu, 1) # data
for(i in 1:10000) {
  ypred[i] ~ dnorm(mu, 1) # predictive depth 1
  zpred[i] ~ dnorm(ypred[i], 1) # predictive depth 2
}

I think this will give 10000 posterior predictive branches, each with a (ypred[i], zpred[i]) pair.

[danielturek]

@perrydv I've made some efficiency improvements to identifying posterior predictive branches. Notably, I've combined the 3 loops into a single loop, which now only contains a single call to getDependencies for each candidate predictive branch node. I've also added testing for the sampler assignments to predictive branches.

Also, the "possibly challenging benchmark case" you presented, the MCMC configuration time is nearly identical for using / not using the option MCMCjointlySamplePredictiveBranches, each are about 22 or 23 seconds on my machine, so I don't see any much of a performance hit.

Thanks for the careful code review. I'd like to say this is ready to merge, but I welcome any further feedback.

[danielturek]

Few more small changes, but again I consider this ready to merge.

[perrydv]

👍 Merging this is good with me.

Here is an idea for where we could go with this and also to avoid conditioning on any posterior predictive nodes in MCMC.

We could give the maps a predictive_IDs similar to (and not mutually exclusive of) top_IDs, latend_IDs and end_IDs. Then getNodeNames could get a predictiveOnly and/or includePredictive option. getDependencies already has an omit option, so one could omit predictive nodes using that. Or it could get a new omitPredictive option to avoid what could become a lot of redundant steps otherwise. Then the calls to getDependencies for every sampler could exclude predictive nodes.

This is still different from grouping predictive nodes into shared branches as this PR does.

[danielturek]

@perry I fully support the idea of adding a predictiveNode flag to the maps. The determination of that could use (via refactoring) the detection logic which this PR adds to the MCMC. For now, I'm merging this PR.

[paciorek]

I also like the predictiveIDs idea. Should we start a NCT issue?