How to calculate credibility interval limits

#1
Hi everyone.
My name is Choi. I'm working on detecting drug adverse reaction signals using
health insurance data. I have a problem in getting credibility interval limits due to
my poor math and statistics below is what I quoted from a paper.
---------------------------------------------------------------------------------
'Ωcan be viewed as the logarithm of the posterior mean of an unknown rate of
incidence μ under the natural assumption that n111 is Po(μ·E111)-distributed with log2μ=Ω and
a gamma prior distribution (or random effects model in a likelihood-based analysis) for μ: G(α,α),
with expected value 1. The choice of prior is made mainly for mathematical convenience, since
due to conjugacy the posterior distribution for μ will also be gamma (but with parameters n111+α
and E111+α, expected value (n111+α)/(E111+α) and variance (n111+α)/(E111+α)^2).
With the Bayesian approach, exact credibility interval limits for μ can be found numerically as
solutions to the following equation, for appropriate posterior quantiles μq :

**(equation 20) **



Specifically, the logarithm of the solutions to (20) for q=0.025 and 0.975, respectively, provides
the upper and lower limits of a two-sided 95 per cent credibility interval for Ω: Ω025 and Ω975.'
----------------------------------------------------------------------------------
I used calculation tools to figure out uq when q=0.025. But they made errors and failed
to get the value. I worder how I can solve this problem.. seeing the table or calculus
tools or whatever.
Thankyou.
 

BGM

TS Contributor
#2
So I just see you are trying to calculate the quantile of the posterior (gamma) distribution. What is your error? Can you provide more information?
 
#3
Thanks for the reply. I wanted to the practical tool to get the value.
Anyway I've got the answer. Below is from the author of the paper.
------------------------------------------------------------------------------------------------------------------------------------------------
For details of how to calculate the credibility intervals from the observed and expected values I would recommend you to read the newly published article: Norén GN, Hopstadius J, Bate A. Shrinkage observed-to-expected ratios for robust and transparent large-scale pattern discovery. Statistical Methods in Medical Research. June 24, 2011



The key in the calculations is that the ratio (O+1/2)/(E+1/2) is Gamma G(O+1/2, E+1/2) distributed making the lower limit of the 95% credibility interval to be calculated from the quantiles of the gamma distribution. As an example, in R (open source statistical software) IC and IC025 may be calculated from the variables observed and expected:

omega <- log2( (observed+0.5)/(expected+0.5) )

omega_025 <- log2(qgamma(p=0.025 , shape=(observed +0.5), rate=(expected +0.5)))

Some implementations (e.g. the quantile function in SAS) uses the parameter scale = 1/rate.



In the referred paper there is an approximation possible to use if you don’t have access to any statistical packages.