***************************************************************************** * Stata do-file to replicate results given in the following presentation: * * TITLE: Running MLwiN from within Stata: The runmlwin command * AUTHORS: George Leckie and Chris Charlton * MEETING: Modern Modeling Methods Conference, University of Connecticut * LOCATION: Bristol * DATE: 26-05-2011 * * George Leckie and Chris Charlton * Centre for Multilevel Modelling, 2011 ***************************************************************************** * The examples in this do-file are adapted from: * * Rabe-Hesketh, S. and Skrondal, A. (2008) Multilevel and Longitudinal * Modeling Using Stata, 2nd Edition, Stata Press. * ***************************************************************************** * To run this do-file in Stata you must be working with: * * (1) Stata 9.2 or higher * * (2) The latest version of MLwiN * * (3) The latest version of the runmlwin command * * runmlwin will inform you if you are working with out of date versions of * Stata or the MLwiN software package. * * You must also tell runmlwin where MLwiN is installed. One way you can do * this is by specifying the MLwiN path by defining a global macro called * MLwiN_path. For example: * * . global MLwiN_path C:\Program Files\MLwiN v2.25\mlwin.exe * * See the "Remarks on installation instructions" section of the runmlwin help * file for further information: * * . help runmlwin * ***************************************************************************** * Change this global macro to point to your MLwiN path global MLwiN_path C:\Program Files (x86)\MLwiN v2.25\mlwin.exe ***************************************************************************** * (1) TWO-LEVEL MULTILEVEL MODELS ***************************************************************************** use "http://www.bristol.ac.uk/cmm/media/runmlwin/tutorial.dta", clear runmlwin normexam cons, /// level2(school: cons) /// level1(student: cons) runmlwin normexam cons, /// level2(school: cons, residuals(u)) /// level1(student: cons) nopause egen pickone = tag(school) preserve keep if pickone==1 egen u0rank = rank(u0) serrbar u0 u0se u0rank, scale(1.96) yline(0) summarize u0 generate u0std = (u0 - r(mean))/r(sd) generate u0uniform = (u0rank - 0.5)/_N generate u0nscore = invnorm(u0uniform) scatter u0std u0nscore, /// yline(0) xline(0) ylabel(-3(1)3) xlabel(-3(1)3) /// aspectratio(1) restore drop u0 u0se generate boy = 1 - girl runmlwin normexam cons standlrt girl, /// level2(school: cons standlrt, residuals(u)) /// level1(student: girl boy, diagonal) nopause test [RP1]var(girl) = [RP1]var(boy) ***************************************************************************** * (2) GROWTH CURVE MODELS ***************************************************************************** use "http://www.stata-press.com/data/mlmus2/asian.dta", clear graph twoway (connect weight age, connect(ascending)), /// ytitle("Weight in Kg") xtitle("Age in years") generate cons = 1 generate age2 = age^2 runmlwin weight cons age age2, /// level2(id: cons age, residuals(u)) /// level1(occ: cons) nopause generate prediction = _b[cons]*cons + _b[age]*age + _b[age2]*age2 /// + u0 + u1*age line prediction age, connect(a) /// ytitle("Predicted weight in Kg") xtitle("Age in years") label define genderlabel 1 "Boy" 2 "Girl" label values gender genderlabel graph twoway (connect weight age, connect(ascending)), by(gender) /// ytitle("Weight in Kg") xtitle("Age in years") generate boy = (gender==1) generate girl = (gender==2) generate boyXage = boy*age generate girlXage = girl*age matrix a = (1,1,1,0,0,1,0,0,1,1) runmlwin weight boy boyXage girl girlXage age2, /// level2(id: boy boyXage girl girlXage, elements(a)) /// level1(occ: boy girl, diagonal) nopause ***************************************************************************** * (3) MULTILEVEL MODELS FOR BINARY RESPONSES ***************************************************************************** use "http://www.stata-press.com/data/mlmus2/guatemala.dta", clear generate cons = 1 runmlwin immun cons kid2p rural pcInd81 , /// level3(cluster: cons) /// level2(mom: cons) /// level1(kid:) /// discrete(distribution(binomial) link(logit) denominator(cons)) /// nopause runmlwin immun cons kid2p rural pcInd81 , /// level3(cluster: cons) /// level2(mom: cons) /// level1(kid:) /// discrete(distribution(binomial) link(logit) denominator(cons) pql2) /// initsprevious maxiterations(40) nopause ***************************************************************************** * (4) SIMULATION STUDIES ARE EASY ***************************************************************************** ***************************************************************************** * (5) MCMC ESTIMATION ***************************************************************************** runmlwin immun cons kid2p rural pcInd81 , /// level3(cluster: cons) /// level2(mom: cons) /// level1(kid:) /// discrete(distribution(binomial) link(logit) denominator(cons)) /// mcmc(on) initsprevious nopause mcmcsum, trajectories mcmcsum [RP2]var(cons), fiveway mcmcsum [RP2]var(cons) runmlwin immun cons kid2p rural pcInd81 , /// level3(cluster: cons) /// level2(mom: cons) /// level1(kid:) /// discrete(distribution(binomial) link(logit) denominator(cons)) /// mcmc(burnin(500) chain(50000) thinning(1) hc(3)) /// initsprevious nopause ***************************************************************************** * (6) CROSS-CLASSIFIED MODELS ***************************************************************************** use "http://www.stata-press.com/data/mlmus2/neighborhood.dta", clear table neighid schid if inrange(neighid,26,38) | inrange(neighid,251,263) /// | inrange(neighid,793,800) gen studentid = _n gen cons = 1 matrix b = (0,.33,.33,.33) runmlwin attain cons, /// level3(schid: cons) /// level2(neighid: cons) /// level1(studentid: cons) /// mcmc(cc) initsb(b) /// nopause matrix b = (0,0,0,0,0,0,0,0,0,.1,.1,.1) gen female = 1 - male runmlwin attain cons p7vrq p7read /// dadocc daded momed dadunemp male deprive, /// level3(schid: cons) /// level2(neighid: cons) /// level1(studentid: cons) /// mcmc(cc) initsb(b) /// nopause ***************************************************************************** * (7) SPATIAL MULTILEVEL MODELS ***************************************************************************** use "http://www.bristol.ac.uk/cmm/media/runmlwin/lips1.dta", clear copy "http://www.bristol.ac.uk/cmm/media/runmlwin/scotdb.dta" "scotdb.dta", replace copy "http://www.bristol.ac.uk/cmm/media/runmlwin/scotcoord.dta" "scotcoord.dta", replace merge area using "scotdb.dta", unique sort drop _merge spmap obs using "scotcoord.dta", id(area) fcolor(Blues) legend(position(10)) clmethod(custom) clbreaks(0 5 10 15 20 40) quietly runmlwin obs cons perc_aff, /// level2(area: cons) /// level1(area:) /// discrete(distribution(poisson) link(log) offset(offs)) nopause runmlwin obs cons perc_aff, /// level2(area:cons) /// level1(area:) /// discrete(distribution(poisson) link(log) offset(offs)) /// mcmc(chain(50000) refresh(500)) initsprevious nopause quietly runmlwin obs perc_aff, /// level2(area: cons) /// level1(area:) /// discrete(distribution(poisson) link(log) offset(offs)) nopause runmlwin obs perc_aff, /// level2(area: cons, carids(neigh1-neigh11) carweights(wcar1-wcar11)) /// level1(cons:) /// discrete(distribution(poisson) link(log) offset(offs)) /// mcmc(chain(50000) refresh(500)) initsprevious nopause quietly runmlwin obs cons perc_aff, /// level3(area: cons) /// level2(area: cons) /// level1(area:) /// discrete(distribution(poisson) link(log) offset(offs)) nopause runmlwin obs cons perc_aff, /// level3(area: cons, carids(neigh1-neigh11) carweights(wcar1-wcar11) residuals(v)) /// level2(area: cons, residuals(u)) /// level1(cons:) /// discrete(distribution(poisson) link(log) offset(offs)) /// mcmc(chain(50000) refresh(500)) initsprevious nopause ***************************************************************************** * (8) EXPORT MODELS TO WINBUGS ***************************************************************************** use "http://www.bristol.ac.uk/cmm/media/runmlwin/tutorial.dta", clear quietly runmlwin normexam cons standlrt, /// level2(school: cons) /// level1(student: cons) /// nopause runmlwin normexam cons standlrt, /// level2(school: cons) /// level1(student: cons) /// mcmc(savewinbugs(model("m.txt", replace) inits("i.txt", replace) data("d.txt", replace))) /// initsprevious nopause view m.txt copy "http://www.bristol.ac.uk/cmm/media/runmlwin/7_3a_WinBUGS_model.txt" "m.txt", replace copy "http://www.bristol.ac.uk/cmm/media/runmlwin/7_3a_WinBUGS_inits.txt" "i.txt", replace view m.txt wbscript , /// model("`c(pwd)'\m.txt") /// inits("`c(pwd)'\i.txt") /// data("`c(pwd)'\d.txt") /// coda("`c(pwd)'\out") /// set(df) /// burn(500) update(5000) /// saving("`c(pwd)'\script.txt", replace) /// quit wbrun, script("`c(pwd)'\script.txt") /// winbugs("C:\WinBUGS14\winbugs14.exe") wbcoda, root("`c(pwd)'\out") clear mcmcsum df, trajectories variables mcmcsum df, variables ***************************************************************************** * (9) WORK EFFICIENTLY ***************************************************************************** ***************************************************************************** * (10) RESOURCES TO HELP YOU LEARN RUNMLWIN ***************************************************************************** ***************************************************************************** exit