95%Credible intervals are smaller than OR

Welcome to the forum for runmlwin users. Feel free to post your question about runmlwin here. The Centre for Multilevel Modelling take no responsibility for the accuracy of these posts, we are unable to monitor them closely. Do go ahead and post your question and thank you in advance if you find the time to post any answers!

Go to runmlwin: Running MLwiN from within Stata >> http://www.bristol.ac.uk/cmm/software/runmlwin/
Post Reply
ManuelDewez
Posts: 16
Joined: Tue Mar 09, 2021 5:54 pm

95%Credible intervals are smaller than OR

Post by ManuelDewez »

Hello,

I am running the following multilevel model with a binary outcome with mcmc:
runmlwin unique_dipstick_avai2cats cons private secondaryhosp paed_mother_child centered_unique_turn19_time centered_exp_capita gov VHIOOP, level2(country2: cons) level1(id:) discrete(distribution(binomial) link(logit) denom(cons) pql2) nopause
runmlwin unique_dipstick_avai2cats cons private secondaryhosp paed_mother_child centered_unique_turn19_time centered_exp_capita gov VHIOOP, level2(country2: cons, residuals(u)) level1(id:) discrete(distribution(binomial) link(logit) denom(cons) pql2) mcmc(burnin(1000) chain(550000)) initsprevious nopause nogroup
I have 504 observations (level 1) and 29 second level units.

My outcome is whether urine dipstick are available (0=no; 1=yes). I am using mcmc because the response proportion is extreme (0=42/504; 1=462/504).

I checked the diagnostics and following the Raftery-Lewis and Brooks-Drappers diagnostics I re run the model with 550,00 iterations.

Some of the odds ratio and Standard deviations are huge for 2 of the variables (OR for gov:1156.948 , and VHIOOP:13730.1), and the 95%Credible interval is smaller than the parameter for VHIOOP: 95%CI = 0.050643 - 164.6105

Burnin = 1000
Chain = 550000
Thinning = 1
Run time (seconds) = 156
Deviance (dbar) = 214.32
Deviance (thetabar) = 195.41
Effective no. of pars (pd) = 18.91
Bayesian DIC = 233.23

Disptick_available OR Std. Dev ESS P [95% Cred.Int]

gov 1156.948 24956.07 8874 0.003 2.891191 2926.618
VHIOOP 13730.1 1120034 14575 0.307 0.050643 164.6105

Can anyone explain why is this happening, and what should I do to get correct OR and 95%CI?

Thank you!

Manuel
ChrisCharlton
Posts: 1348
Joined: Mon Oct 19, 2009 10:34 am

Re: 95%Credible intervals are smaller than OR

Post by ChrisCharlton »

That does sound a bit strange. Have you tried plotting the chains? What do the results look like it you don't convert them to odds-ratios?
ManuelDewez
Posts: 16
Joined: Tue Mar 09, 2021 5:54 pm

Re: 95%Credible intervals are smaller than OR

Post by ManuelDewez »

Hello Chris,

The ACF plot shows quite some autocorrelation and the MCSE plot show that I should run the model with more iterations ( which I did)
The non-ORs results look fine ( the parameters are within the 95% Credible Intervals).

Gov and VHIOOP are 2 categories of the same variable. I think the issue is that there are too few observations for this variable (only 8 for gov and 5 for VHIOOP:


unique_dip |
stick_avai | source_financing2
2cats | compshi gov vhi_oop | Total
-----------+---------------------------------+----------
0 | 29 8 5 | 42
1 | 221 211 30 | 462
-----------+---------------------------------+----------
Total | 250 219 35 | 504

And this is creating the troubles ?
I think I will just ignore this variable for this particular model (I don't have this issue with other models where there is a more balanced distribution of the observations across the outcome's binary categories). What do you think?

Thank you.

Manuel
ChrisCharlton
Posts: 1348
Joined: Mon Oct 19, 2009 10:34 am

Re: 95%Credible intervals are smaller than OR

Post by ChrisCharlton »

The way that the credible intervals for this would be calculated is as follows:
  1. Fetch the MCMC chain for the parameter
  2. Exponentiate the values to convert to odds-ratios
  3. Sort the chain
  4. Extract the 2.5% and 97.5% quantiles
I suppose that if the chain contained a few extreme values then this could result in the mean being outside the range of the credible intervals. You might also want to look at the mode and median values for comparison. For a more visual examination of the chain distribution you can use the fiveway option in the mcmcsum command to look at the kernel density.
ManuelDewez
Posts: 16
Joined: Tue Mar 09, 2021 5:54 pm

Re: 95%Credible intervals are smaller than OR

Post by ManuelDewez »

Thank you Chris,

this is the syntax I am using :
runmlwin unique_dipstick_avai2cats cons private secondaryhosp paed_mother_child centered_unique_turn19_time centered_exp_capita gov VHIOOP if unique_hospital==1, level2(country2: cons) level1(id:) discrete(distribution(binomial) link(logit) denom(cons) pql2) nopause
followed by :
runmlwin unique_dipstick_avai2cats cons private secondaryhosp paed_mother_child centered_unique_turn19_time centered_exp_capita gov VHIOOP if unique_hospital==1, level2(country2: cons, residuals(u)) level1(id:) discrete(distribution(binomial) link(logit) denom(cons) pql2) mcmc(burnin(1000) chain(550000)) initsprevious nopause nogroup
to get OR I use:
runmlwin, or

this provides ORs and 95 %Cred. Int.


How do I sort and extract the 2.5 and 97.5 % quantiles?
Thank you
ChrisCharlton
Posts: 1348
Joined: Mon Oct 19, 2009 10:34 am

Re: 95%Credible intervals are smaller than OR

Post by ChrisCharlton »

You can obtain the MCMC chains as a dataset with the following command:

Code: Select all

mcmcsum, getchains
Once you have this you can use the standard Stata data manipulation commands to process it, i.e. gen, sort and centile .

You ought to get the same results as before as this should be equivalent to what runmlwin is doing itself.
ManuelDewez
Posts: 16
Joined: Tue Mar 09, 2021 5:54 pm

Re: 95%Credible intervals are smaller than OR

Post by ManuelDewez »

That worked well, thank you Chris!
Post Reply