From treatment responses to a predictive decision

Question

Is the response rate plausibly above 50%?

Seven of ten participants responded. We use beta(2,2), a symmetric prior that places little weight at implausible extremes. The model is intentionally conjugate so its beta(9,5) posterior can validate the MCMC workflow.

data { int n; int y; }
parameters { real<lower=0, upper=1> theta; }
model {
  theta ~ beta(2, 2);
  y ~ binomial(n, theta);
}
generated quantities {
  int y_rep = binomial_rng(n, theta);
}

Workflow

Compile once, sample four streams

Map<String, double[]> data = new LinkedHashMap<>();
data.put("n", new double[] {10});
data.put("y", new double[] {7});
CompiledModelScript script = ModelScript.compile(source, data);

SamplingOptions options = SamplingOptions.builder()
    .warmupIterations(500).sampleIterations(1000)
    .targetAcceptance(0.85).build();
ChainResult[] chains = Chains.parallel(new NoUTurnSampler(),
    script.model(), new double[][] {{-1}, {-.25}, {.25}, {1}},
    options, 41024L, 4);

Before sampling, the production workflow also calls Gradients.check. After sampling, it rejects any failed chain and reviews diagnostics and plots before computing decisions.

Posterior scale

Transform each draw before answering the question

int aboveHalf = 0;
double posteriorMean = 0.0;
for (ChainResult chain : chains) {
  for (int draw = 0; draw < chain.size(); draw++) {
    double theta = script.model().state(chain.sample(draw))
        .scalar("theta");
    posteriorMean += theta;
    if (theta > 0.5) aboveHalf++;
  }
}
int total = chains.length * chains[0].size();
posteriorMean /= total;
double probabilityAboveHalf = aboveHalf / (double) total;

Prediction

Simulate a future cohort

Map<String, double[]> generated = script.generate(
    chains[0].sample(0), new MersenneTwister(991));
double futureResponses = generated.get("y_rep")[0];

Repeat generation over posterior draws to obtain the posterior-predictive distribution for ten future participants. This includes uncertainty in the response probability; plugging only the posterior mean into a binomial would understate that uncertainty.