3
$\begingroup$

If one is calculating odds ratio with a,b,c and d counts, I believe variance of log(OR) is given by

var_log_OR = (1/a + 1/b + 1/c + 1/d)

Hence one can calculate 95% confidence intervals of OR as follows:

SE_log_OR = sqrt(var_log_OR)

CI_lower_log_OR = log(OR) - 1.96*SE_log_OR
CI_upper_log_OR = log(OR) + 1.96*SE_log_OR

CI_lower_OR = exp(CI_lower_log_OR)
CI_upper_OR = exp(CI_upper_log_OR)

But how can we calculate SE of OR?

$\endgroup$
5
$\begingroup$

@FrankHarrell is right that the standard error for an odds ratio is a problematic number in the sense that you can do better by testing on the corresponding log(odds ratio) scale, as the sampling distribution of the log(odds ratio) is more likely to be normally distributed.

Nonetheless, the standard error of the odds ratio does exist, even if it is not that useful. One possible estimate is to use the delta method to move from the standard error of the log(odds ratio) to an approximation of the standard error of the odds ratio.

$\sqrt{(1/a + 1/b + 1/c + 1/d)}\times\frac{a\times d}{b\times c}$

| cite | improve this answer | |
$\endgroup$
5
$\begingroup$

OR is not a valid quantity to compute a SE of in the sense that it cannot have a symmetric distribution. Applying +/- SE to it may lead to negative ORs.

| cite | improve this answer | |
$\endgroup$
  • $\begingroup$ But it is reported in output of meta-analysis package metafor as can be seen in output of rma() function in answer to this question: stats.stackexchange.com/questions/156480/… $\endgroup$ – rnso Jun 12 '15 at 3:10
  • 2
    $\begingroup$ I hope that what is stated there is the log OR otherwise it needs to be fixed. $\endgroup$ – Frank Harrell Jun 12 '15 at 3:14
  • $\begingroup$ Can we simply say that (CI_upper_OR - CI_lower_OR)/(2*1.96) will approximate SE of OR ? $\endgroup$ – rnso Jun 12 '15 at 3:37
  • $\begingroup$ No, since upper - OR is not equal to OR - lower so you are averaging two unequal quantities. $\endgroup$ – Frank Harrell Jun 12 '15 at 11:17
  • $\begingroup$ @rnso No, metafor is computing everything on the log-OR scale. The answer given to the question you linked to is misleading. $\endgroup$ – Wolfgang Jun 14 '15 at 7:24
2
$\begingroup$

Instead of a standard error why not compute the standard deviation of the posterior distribution of the OR? You can solve for it numerically very easily using an MCMC sampler.

Here is some R and JAGS code to do so.

################################################################
###                                                          ###
###      Contingency Table Analysis for Obestity Data        ###
###                                                          ###
################################################################

# Required Pacakges
library("ggplot2")
library("runjags")
library("parallel") # sets parallelization for MCMC


# set up the model
mod = 'model {

################################################
###    Greater than 80th BMI Percentile      ###
################################################

  # marginal likelihood functions
  n11_G ~ dbin(pi_one_G, n1_plus_G)
  n21_G ~ dbin(pi_two_G, n2_plus_G)

  #priors
  pi_one_G ~ dbeta(1,1)
  pi_two_G ~ dbeta(1,1)

  # transformations
  rho_G <- pi_two_G/pi_one_G
  theta_G <- pi_two_G*(1-pi_one_G)/(pi_one_G*(1-pi_two_G))
  delta_G <- pi_two_G-pi_one_G



}'

# set up the data
Dat = list(n1_plus_G = 108,
           n2_plus_G = 88,
           n11_G = 68,
           n21_G = 44)

# Monitor these variables
Vars = c("pi_one_G","pi_two_G","rho_G","theta_G","delta_G")


# set up MCMC parameters
inits1=list(.RNG.name= "base::Wichmann-Hill",
            .RNG.seed= 12341)
inits2=list(.RNG.name= "base::Marsaglia-Multicarry",
            .RNG.seed= 12342)
inits3=list(.RNG.name= "base::Super-Duper",
            .RNG.seed= 12343)
inits4=list(.RNG.name= "base::Mersenne-Twister",
            .RNG.seed= 12344)

chains = 4
burn = 5000
samp = 10000
adapt = 5000
thin = 1

# parallel chains
cl = makeCluster(4)

# MCMC estimation
HjagsOut = run.jags(model = mod, monitor = Vars, data=Dat, n.chains=chains, thin = thin,
                    burnin = burn, sample = samp, adapt=adapt, method="rjparallel",method.options=list(cl=cl),
                    inits=list(inits1,inits2,inits3,inits4))

#summarize results
summary(HjagsOut)

plot(HjagsOut, layout=c(4,2))

The odds ratio parameter ($\theta_\text{G}$) is simply a function of the samples from the binomial parameters. This of course assumes certain study design. In this case I was looking at the difference in children's BMI percentile group (80th and above or below 80th) from a control and experimental group, pre and post intervention treatment. Therefore, the row totals (number of children in the experimental group and control group respectively) were fixed.

The beta(1,1) prior is equivalent to a Uniform(0,1) prior and could easily be changed to the Jeffreys's beta(0.5,0.5) prior or anything you desire.

| cite | improve this answer | |
$\endgroup$
  • $\begingroup$ That is interesting suggestion. Can you show briefly how to compute it (even if in R code). $\endgroup$ – rnso Jun 12 '15 at 3:35
  • $\begingroup$ Sure I will dig up some code. $\endgroup$ – Nathan L Jun 12 '15 at 3:49
  • $\begingroup$ Of course you can do that but the usefulness of it is questionable. How do you interpret the SD of an asymmetric distribution? Worse is the dichotomization of BMI. $\endgroup$ – Frank Harrell Jun 12 '15 at 11:19
  • $\begingroup$ Well he was interested in the SD. I would suggest finding the (1-alpha) HPD interval and perhaps probabilites that the OR is within a specific interval of interest. As for the dichotomization of BMI, I agree, but the PIs on this project were interested in looking at their data this way initially. In the actual paper (still under peer review) I expand beyond the 2x2 table view of the data and examine the continuous BMI response as well. $\endgroup$ – Nathan L Jun 12 '15 at 14:07
  • $\begingroup$ The initial look using dichotomized BMI will be nothing but misleading. HPD interval is a good idea, and avoids SD. $\endgroup$ – Frank Harrell Jun 12 '15 at 18:12

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.