What is the most accurate way to characterize probability from a posterior?

I have a philosophical question about interpreting samples from a posterior. Consider the simple example:

``````library(rstanarm)
set.seed(1839)
x     <- rnorm(100)
y     <- x + rnorm(100, 0, 4)
dat   <- data.frame(x, y)
model <- stan_glm(
y ~ x,
prior = normal(0, 2),
prior_intercept = normal(0, 2),
prior_aux = cauchy(0, 5)
)
``````

Let’s say my client wants to know if `x` positively predicts `y`. I could calculate how many of my posterior samples were above zero:

``````mean(as.matrix(model)[, "x"] > 0)
``````

Which returns `0.997`.

Which would you all consider to be an accurate characterization of the probability? Consider the interpretations:

A) “The probability that there is a positive relationship between x and y is .997”

B) “The probability that there is a positive relationship between x and y, given the data, is .997”

C) “The posterior probability that there is a positive relationship between x and y, is .997” (Problem here is that clients won’t know what this means, likely)

D) “The probability that there is a positive relationship between x and y, given our model and priors, is .997”

E) “Given the model and priors we specified, as well as the data we observed, there is a .997 probability that there is a positive relationship between x and y”

(Nevermind for now that an estimate of .0001 might as well be zero, but is still considered “positive” for the time being).

I really enjoy MCMC and sampling from the posterior, because it leads to intuitive probability statements. However, I do not want to overstate the results, either. What would you consider to be the proper way to interpret proportions form a posterior in probabilistic terms?

Other ways you could look at the data, which I assume would follow the same probability statements:

``````> mean(as.matrix(model)[, "x"] > 0.5)
[1] 0.9415
> mean(as.matrix(model)[, "x"] > 0.5 & as.matrix(model)[, "x"] < 1.5)
[1] 0.7185
> mean(as.matrix(model)[, "x"] < 0.5 & as.matrix(model)[, "x"] < 1.5)
[1] 0.0585``````

I would go with

Under our model (which consists of a definition of the parameter space, a prior distribution over that parameter space, and one or more distributions for how data are generated by the parameters), given the data, the probability that `x` is greater than zero is essentially 1.

The reason I would say essentially 1 or greater than 0.99 is that you are not going to get three decimal places of precision from 4000 draws, but probably no one cares about that level of precision anyway.

In terms of your five statements,

(A) I think Bayesians would take this to mean a posterior probability, but in principle you could say “the probability that there is a positive relationship between x and y is 0.997” under the prior. Also, any Bayesian should understand that the statement is with respect to a particular model rather than the result of integrating over models unless you specifically refer to integrating over models.

(B) and © are the same. (D) could be interpreted as a prior probability and (E) is pretty close to how I would phrase it, although I tend to think of the declaration of the parameter space and the priors to be part of the model. Non-Bayesians have a tendency to think of “the model” as “the likelihood” or “the procedure for minimizing a loss function”.

In the end, I doubt the semantics matter much. People who are not versed in this stuff are liable to misinterpret whatever you say anyway. But such people are more likely to interpret results in a roughly Bayesian way than in a frequentist way, which is a problem when they are listening to frequentists.

3 Likes

Thanks for your answer! I think that “Under the model, the probability that x positively predicts y is greater than .99” is a reasonable statement, based on what you’ve described. \

This is because the colloquial definition of what a “model” is includes all of the “statistical things” we do (parameters, priors, likelihoods, etc.). Only within more statistical circles would there be some confusion as to what a “model” means (e.g., only the likelihood).