re: Stata Command to estimate Bayesian Interval for the risk of rare complications

[Edsel Ing]

Is there a Stata function / plug-in to estimate the risk of rare medical complications, i.e. Bayesian central interval?

Estimating the risk of rare complications: is the ‘rule of three’ good enough? Ludbrook, John ; Lew, Michael J ANZ journal of surgery, 2009-07, Vol.79 (7-8), p.565-570

p for the sample binomial proportion, x for the numerator of that
proportion and n (sample size) for the denominator. Thus, p = x/n.

Thanks

[Leonardo Guizzetti]

Here is an implementation using -bayesmh- which uses the same model specified in that paper (though there is no clear discussion of the type of Bayesian credible interval being used -- presumably it is the equal-tailed interval).

clear *
cls
set seed 17

* e = events, n=total observations
local n 100
local e 5

* convert this into a dataset suitable for -bayesmh-
set obs 2
gen byte y = _n-1
gen int n = `e' if y==1
replace n = `n' - `e' if y==0
list

* run the Bayesian model
bayesmh y [fw=n], likelihood(dbern({pr})) prior({pr}, beta(1,1))

On the practicality of the approach, I have several concerns with the paper and suggest that you take it with a grain of salt. When working with small N and rare outcomes, there are no good methods because there is relatively little information from which to make an inference. But if N is not so small (N>50 say), the Wilson (score) interval is reasonable to use. The comments about the Clopper-Pearson interval are well known, especially its sawtooth pattern and conservative coverage.

The authors make much of the Bayesian approach without really getting into the fact that a Bayesian interpretation of a credible interval does not have the same interpretation as the confidence interval, yet are happy to accept them as the same. I don't think there's anything wrong with picking frequentist or Bayesian as a philosophical stance, but comparing them in this way seems unfair and is a bit like comparing apples to oranges.

The Bayesian model of a binomial parameter with a uniform prior is *not* itself "uninformative", but the authors skip over its justification entirely, though this is the common refrain for its use. A uniform (flat) prior in the context of a binomial parameter is mathematically identical to a beta(1,1) distribution because it must be constrained to [0, 1]. This has the effect of assuming a prior with 2 observations -- 1 event and 1 non-event -- so this prior adds information to the existing data. When N is large or the probability parameter is not so extreme, this doesn't matter because the observed data overwhelm the prior, but this is not always true if either assumption is violated. This is why the Bayesian interval never contains 0 even when the true parameter value is 0. The authors might have commented on this if they examined bias, as well as coverage, which can be just as important for examining safety signals with rare events.

[Nick Cox]

I'll add that the Jeffreys interval available under ci proportions, jeffreys was originally proposed in Bayesian spirit but also has a frequentist interpretation. I find that it works well across a range of sample sizes.

I've not updated it since 2003 but there are some relevant references in

ssc type cij.hlp

that didn't make it to the manual entry.