Hi, I am new to runmlwin and am trying to obtain the DIC for a multinomial random effects model in Stata 11 for Windows. I came across runmlwin while searching for a way to compute DIC in Stata following a gllamm with mlogit link. Other than runmlwin, I do not know of any other way to compute DIC in Stata.
First I decided to run the model using example data. The model runs successfully however it does not produce a DIC estimate.
I tried estout and while the command seems to be running, it does not produce a DIC estimate. I would appreciate any help on this issue. Here is the code I used:
use http://www.bristol.ac.uk/cmm/media/runmlwin/bang, clear
generate lc1 = (lc==1)
generate lc2 = (lc==2)
generate lc3plus = (lc>=3)
runmlwin use4 cons lc1 lc2 lc3plus, level2(district: cons) level1(woman) discrete(distribution(multinomial) link(mlogit) denom(cons) basecategory(4)) nopause
estout, stats(dic)
=====================================================================================================================================
Bogdan Cristescu
Ph.D. Ecology Candidate
Department of Biological Sciences
University of Alberta
Canada
DIC statistic for multinomial random effects model
-
bogdan2012
- Posts: 3
- Joined: Mon Jun 04, 2012 8:16 pm
-
GeorgeLeckie
- Site Admin
- Posts: 432
- Joined: Fri Apr 01, 2011 2:14 pm
Re: DIC statistic for multinomial random effects model
Hi Bogdan,
The DIC statistic is essentially the Bayesian equivalent of the AIC statistic which you might be used to reporting after fitting models by maximum likelihood.
The DIC is therefore only relevant after fitting models by MCMC.
The following syntax fits your model by MCMC where we use the IGLS quasilikelihood estimates as starting values for the MCMC chains
Here is the associated output and at the very end you can see that the DIC statistic for this model is 5835.916
Best wishes
George
The DIC statistic is essentially the Bayesian equivalent of the AIC statistic which you might be used to reporting after fitting models by maximum likelihood.
The DIC is therefore only relevant after fitting models by MCMC.
The following syntax fits your model by MCMC where we use the IGLS quasilikelihood estimates as starting values for the MCMC chains
Code: Select all
use http://www.bristol.ac.uk/cmm/media/runmlwin/bang, clear
generate lc1 = (lc==1)
generate lc2 = (lc==2)
generate lc3plus = (lc>=3)
runmlwin use4 cons lc1 lc2 lc3plus, ///
level2(district: cons) ///
level1(woman) ///
discrete(distribution(multinomial) link(mlogit) denom(cons) basecategory(4)) ///
nopause
runmlwin use4 cons lc1 lc2 lc3plus, ///
level2(district: cons) ///
level1(woman) ///
discrete(distribution(multinomial) link(mlogit) denom(cons) basecategory(4)) ///
mcmc(on) initsprevious ///
nopause
estout, stats(dic)
exitCode: Select all
. use http://www.bristol.ac.uk/cmm/media/runmlwin/bang, clear
. generate lc1 = (lc==1)
. generate lc2 = (lc==2)
. generate lc3plus = (lc>=3)
. runmlwin use4 cons lc1 lc2 lc3plus, level2(district: cons) level1(woman) discrete(distribution(multinomial) link(mlogit) denom(cons)
> basecategory(4)) nopause
MLwiN 2.25 multilevel model Number of obs = 2867
Unordered multinomial logit response model
Estimation algorithm: IGLS, MQL1
-----------------------------------------------------------
| No. of Observations per Group
Level Variable | Groups Minimum Average Maximum
----------------+------------------------------------------
district | 60 3 47.8 173
-----------------------------------------------------------
----------------------------------
Contrast | Log-odds
-------------+--------------------
1 | 1 vs. 4
2 | 2 vs. 4
3 | 3 vs. 4
----------------------------------
Run time (seconds) = 5.37
Number of iterations = 10
------------------------------------------------------------------------------
| Coef. Std. Err. z P>|z| [95% Conf. Interval]
-------------+----------------------------------------------------------------
Contrast 1 |
cons_1 | -3.985474 .3137803 -12.70 0.000 -4.600472 -3.370476
lc1_1 | 2.150926 .3391147 6.34 0.000 1.486273 2.815579
lc2_1 | 2.689987 .3312585 8.12 0.000 2.040732 3.339242
lc3plus_1 | 2.658316 .3145481 8.45 0.000 2.041813 3.274819
-------------+----------------------------------------------------------------
Contrast 2 |
cons_2 | -1.588393 .123751 -12.84 0.000 -1.830941 -1.345846
lc1_2 | .7063671 .1435437 4.92 0.000 .4250267 .9877075
lc2_2 | .6866727 .1518692 4.52 0.000 .3890146 .9843309
lc3plus_2 | .2447693 .1307283 1.87 0.061 -.0114533 .500992
-------------+----------------------------------------------------------------
Contrast 3 |
cons_3 | -2.577769 .1696013 -15.20 0.000 -2.910182 -2.245357
lc1_3 | .7263379 .2172708 3.34 0.001 .3004949 1.152181
lc2_3 | 1.061379 .2126313 4.99 0.000 .6446293 1.478129
lc3plus_3 | 1.125481 .1776694 6.33 0.000 .7772552 1.473707
------------------------------------------------------------------------------
------------------------------------------------------------------------------
Random-effects Parameters | Estimate Std. Err. [95% Conf. Interval]
-----------------------------+------------------------------------------------
Level 2: district |
var(cons_1) | .3494724 .1123334 .1293029 .5696419
cov(cons_1,cons_2) | .1105162 .0695678 -.0258342 .2468666
var(cons_2) | .2887423 .0840008 .1241038 .4533808
cov(cons_1,cons_3) | .027815 .0726465 -.1145695 .1701995
cov(cons_2,cons_3) | -.0408545 .0636466 -.1655994 .0838905
var(cons_3) | .2602587 .0939646 .0760914 .444426
------------------------------------------------------------------------------
. runmlwin use4 cons lc1 lc2 lc3plus, level2(district: cons) level1(woman) discrete(distribution(multinomial) link(mlogit) denom(cons)
> basecategory(4)) mcmc(on) initsprevious nopause
MLwiN 2.25 multilevel model Number of obs = 2867
Unordered multinomial logit response model
Estimation algorithm: MCMC
-----------------------------------------------------------
| No. of Observations per Group
Level Variable | Groups Minimum Average Maximum
----------------+------------------------------------------
district | 60 3 47.8 173
-----------------------------------------------------------
----------------------------------
Contrast | Log-odds
-------------+--------------------
1 | 1 vs. 4
2 | 2 vs. 4
3 | 3 vs. 4
----------------------------------
Burnin = 500
Chain = 5000
Thinning = 1
Run time (seconds) = 115
Deviance (dbar) = 5730.03
Deviance (thetabar) = 5624.14
Effective no. of pars (pd) = 105.89
Bayesian DIC = 5835.92
------------------------------------------------------------------------------
| Mean Std. Dev. ESS P [95% Cred. Interval]
-------------+----------------------------------------------------------------
Contrast 1 |
cons_1 | -4.059774 .3121117 39 0.000 -4.66219 -3.505688
lc1_1 | 2.127718 .3477726 48 0.000 1.449625 2.757305
lc2_1 | 2.687625 .330258 47 0.000 2.084104 3.309828
lc3plus_1 | 2.634509 .3202318 40 0.000 2.042424 3.246418
-------------+----------------------------------------------------------------
Contrast 2 |
cons_2 | -1.713854 .1354397 94 0.000 -2.000681 -1.465363
lc1_2 | .7865393 .1502616 304 0.000 .5004651 1.086808
lc2_2 | .7747856 .1511993 302 0.000 .4740293 1.077601
lc3plus_2 | .3000943 .1258711 247 0.007 .0635232 .5542534
-------------+----------------------------------------------------------------
Contrast 3 |
cons_3 | -2.730865 .1761592 107 0.000 -3.118661 -2.402368
lc1_3 | .8085474 .2213781 203 0.000 .3812687 1.257318
lc2_3 | 1.151097 .2272489 207 0.000 .7075114 1.640513
lc3plus_3 | 1.205638 .1776815 140 0.000 .8746599 1.565538
------------------------------------------------------------------------------
------------------------------------------------------------------------------
Random-effects Parameters | Mean Std. Dev. ESS [95% Cred. Int]
-----------------------------+------------------------------------------------
Level 2: district |
var(cons_1) | .5700164 .1726627 312 .3042242 .9617038
cov(cons_1,cons_2) | .324783 .1126461 470 .1382326 .5717295
var(cons_2) | .4244059 .1223081 474 .2307351 .6926774
cov(cons_1,cons_3) | .2580462 .1155609 208 .0628227 .5205603
cov(cons_2,cons_3) | .1424723 .0915538 471 -.0210082 .3387358
var(cons_3) | .3873758 .134996 333 .1811532 .7002534
------------------------------------------------------------------------------
. estout, stats(dic)
-------------------------
.
b
-------------------------
FP1
cons_1 -4.059774
lc1_1 2.127718
lc2_1 2.687625
lc3plus_1 2.634509
-------------------------
FP2
cons_2 -1.713854
lc1_2 .7865393
lc2_2 .7747856
lc3plus_2 .3000943
-------------------------
FP3
cons_3 -2.730865
lc1_3 .8085474
lc2_3 1.151097
lc3plus_3 1.205638
-------------------------
RP2
var(cons_1) .5700164
cov(cons_~2) .324783
var(cons_2) .4244059
cov(cons_~3) .2580462
cov(cons_~3) .1424723
var(cons_3) .3873758
-------------------------
RP1
cov(P\_P) 1
-------------------------
dic 5835.916
-------------------------
George
-
bogdan2012
- Posts: 3
- Joined: Mon Jun 04, 2012 8:16 pm
Re: DIC statistic for multinomial random effects model
Thank you George that is very helpful.
I have another question that is not directly related to runmlwin but still interesting because of its broad applicability. My understanding is that AIC or AICc are not reliable for ranking candidate maximum likelihood models with mixed effects (such as random intercept, or random intercept and random slope). Some authors use R packages to rank such models based on DIC. I wonder if there is some confusion in the literature over which IC is suitable for ranking mixed models? I would be interested to learn what different people use and to what extent they think their IC of choice is appropriate.
Thanks in advance,
Bogdan
=====================================================================================================================================
Bogdan Cristescu
Ph.D. Ecology Candidate
Department of Biological Sciences
University of Alberta
Canada
I have another question that is not directly related to runmlwin but still interesting because of its broad applicability. My understanding is that AIC or AICc are not reliable for ranking candidate maximum likelihood models with mixed effects (such as random intercept, or random intercept and random slope). Some authors use R packages to rank such models based on DIC. I wonder if there is some confusion in the literature over which IC is suitable for ranking mixed models? I would be interested to learn what different people use and to what extent they think their IC of choice is appropriate.
Thanks in advance,
Bogdan
=====================================================================================================================================
Bogdan Cristescu
Ph.D. Ecology Candidate
Department of Biological Sciences
University of Alberta
Canada
-
GeorgeLeckie
- Site Admin
- Posts: 432
- Joined: Fri Apr 01, 2011 2:14 pm
Re: DIC statistic for multinomial random effects model
Hi Bogdan,
Interesting point, but you're right one that has broader relevance that just runmlwin. Suggest you pose it on the multilevel modelling discussion list
https://www.jiscmail.ac.uk/cgi-bin/weba ... multilevel
as you will likely get some good references to read from the user community.
George
Interesting point, but you're right one that has broader relevance that just runmlwin. Suggest you pose it on the multilevel modelling discussion list
https://www.jiscmail.ac.uk/cgi-bin/weba ... multilevel
as you will likely get some good references to read from the user community.
George
-
bogdan2012
- Posts: 3
- Joined: Mon Jun 04, 2012 8:16 pm
Re: DIC statistic for multinomial random effects model
Thanks it is now posted on the list you suggested.
Bogdan
Bogdan