Again, manipulating the joint density isn’t sufficient without considering the appropriate expectation. In particular because the order statistic conditions on a prior sample one has to be careful to respect that conditioning otherwise the integrals will remain intractable.
In looking back on this again I realized that I was trying to prove the wrong thing. The marginal conditional
corresponds to simulating a single model configuration from the prior but then multiple observations and then multiple posterior samples before computing the rank, which is not how we construct the SBC algorithm in the paper. What we actually do is sample a single observation for each model configuration and then a single observation before generating multiple posteriors samples and then computing the rank. This subtle difference implies that the probabilities in the order statistic are actually conditional on the data as well,
with an extra integral over y on the outside.
This conditional structure of the joint then allows us to simplify \pi(\theta' \mid \theta, y) = \pi(\theta' | y) from which the proof proceeds basically as intended.
I’m attaching what I believe is a proper proof for anyone interested.sbc_calcs_proper.pdf (96.8 KB)