MCMC Macro

Welcome to the forum for MLwiN users. Feel free to post your question about MLwiN software 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!

Remember to check out our extensive software FAQs which may answer your question: http://www.bristol.ac.uk/cmm/software/s ... port-faqs/
yongjookim78
Posts: 45
Joined: Tue Jan 10, 2017 3:36 am

MCMC Macro

Post by yongjookim78 »

Hello,

I have just started to use MLwiN Macro. Mostly I am learning macro commands by checking what appears on the Command Interface window after I run a procedure by manually (clicking button by using mouse). It has worked well except for running models. For IGLS estimation, "STAR" command worked. However, for MCMC, when I used the below commands, the output was different than what I get from when I actually click "START".

BATCh 1
MAXI
STARt

I would appreciate your wisdom!

Best wishes,
YongJoo
ChrisCharlton
Posts: 1250
Joined: Mon Oct 19, 2009 10:34 am

Re: MCMC Macro

Post by ChrisCharlton »

The start command corresponds to estimation via the (R)IGLS estimation method. To run the model with MCMC you will need to use the MCMC command twice, once for the setup and burnin, and then for the main iterations. The MLwiN help file describes the command as follows (when entering the command you would only include the items in bold, replacing the value with settings appropriate for your model.):

MCMC

This command performs MCMC estimation for the currently set up model. It is advisable to run an IGLS or RIGLS run before running the MCMC methods. There is a different syntax for the burn-in and the main run of the simulation. A burn-in (of length 0 if necessary) MUST be run first to set up the model.

The syntax for the burn-in is as follows:

MCMC <0> burn in for <value 1> iterations, use adaptation method <value 2>, scale factor <value 3>, %acceptance <value 4>, tolerance <value 5> {residual starting values in C1, s.e. residual starting values) in C2} {priors in C3} fixed effect method <value 6> residual method <value 7> level 1 variance method <value 8> other levels variance method <value 9> default prior <value 10> model type <value 11>.

This command starts the MCMC sampler and burns in for <value 1> iterations. The settings for <value 2> to <value 5> will only come into effect if some parameters are updated by MH sampling but must be given.

Value 1 = no of burn in iterations.
Value 2 = 1 if adaptation is to be used, 0 otherwise.
Value 3,Value 4 and Value 5 are explained in the help system (See under MH Settings). Value 4 and Value 5 are ignored if Value 2 = 0).

The user has the option of supplying starting values for residuals at level 2 and above. These should be stacked into 1 column C1, higher level residuals first. If you have more than one set of residuals at a level then residuals corresponding to explanatory variables with lower number get stacked first. If you choose to supply residual starting values then you must stack residuals for every random variable at level 2 and above into the input column.

C2 is required if C1 is specified. This column should contain the standard error matrices for the residuals. These should be in the same format as the standard errors output from the RESI command when the RCOV command takes value 2 i.e. if there is more than one parameter random at a level then the complete (lower triangle) covariance matrix for the residuals at that level is required. The stacking should be in the same level order as for C1.

The user also has the option is to supply prior information (C3) for the parameter estimates. The components for which prior information can be specified differ in the fixed and random parts of the model. In the fixed part prior information can be supplied on a per parameter basis. In the random part of the model prior information is displayed on a per level basis, that is prior information must be given for the full covariance matrix at a given level. For fixed parameters prior information takes the form of a mean and SD, for covariance matrices prior information takes the form of estimate values and a sample size. See the section on "Priors" in the MLwiN help system for more details on the meaning of these terms. The format that the prior column should take is described in the PRIOr command.

As the MCMC command is generic <value 6>, <value 7>, <value 8> and <value 9> need to be input for the sampler to know which methods to use for each set of parameters. These parameters take the values 1 for Gibbs Sampling, 2 for univariate MH Sampling and 3 for multivariate MH Sampling. Note that not all combinations of methods are available for all sets of parameters and all models.

<value 10> indicates which default priors are being used for variance parameters. This parameter takes the value 1 for Gamma priors or 0 for Uniform on the variance scale priors. See the section on "Priors" in the MLwiN help system for more details on the meaning of these default priors.

<value 11> indicates the model type being fitted. This parameter takes the values 1 for Normal models, 2 for Binomial models, 3 for Poisson models and 4 for multivariate Normal models.

The syntax for the main chain is as follows:

MCMC <1> continue for <value 1> (stored) iterations, thinning every <value 2> iterations, appending stacked trace to C1, likelihood to C2, mean of all updates to C3, standard deviation of all updates to C4, run <value 3> Metropolis Hastings updates per iteration for model type <value 4>.

The MCMC 1 command continues with the settings used in the last MCMC 0 command. The sampler will run for <value 1>*<value 2> iterations but only store every <value 2>th iteration to store a total of <value 1> iterations. These <value 1> iterations are stored in C1, an iteration at a time in the order fixed parameters followed by random parameters, highest level first. The likelihood column holds the likelihood value for each stored update. The final two columns are the means and standard deviations of parameters taken from all updates irrespective of the thinning setting. For each iteration any steps run by MH sampling will be updated <value 3> times, whilst any Gibbs steps are run just once. Model type <value 4> is identical to the MCMC 0 command.
yongjookim78
Posts: 45
Joined: Tue Jan 10, 2017 3:36 am

Re: MCMC Macro

Post by yongjookim78 »

Hi,

Thank you for your kind help!

Actually I am very beginner of using MLwiN Commands as well as MCMC method. The one that I am using is Metropolis-Hastings algorithm with non-informative priors. So, manually (without using commands interface), I have run models first by running IGLS (automatically to set up starting values), and then running MCMC without designating priors or anything further than the default option.

So, my commands are two parts:

1. For IGLS:

tole 2
MTOT 5000
EMODe 3
EMODe 0
METH 1
EMODe 0
METH 1

[model specification]

START

2. For MCMC:

tole 2
MTOT 5000
EMODe 3
EMODe 3

BATCh 1
MAXI
START

Could you provide a bit more comments on my codes here? I tried to follow the manual that you've kindly attached, but honestly it is a bit hard for me to follow and apply at this moment. I would appreciate your help!!

Best wishes,
YongJoo
ChrisCharlton
Posts: 1250
Joined: Mon Oct 19, 2009 10:34 am

Re: MCMC Macro

Post by ChrisCharlton »

I have added comments in bold to the code that you provided below:

tole 2 - set (R)IGLS convergence tolerance to 2
MTOT 5000 - set desired length of MCMC chain to 5000
EMODe 3 - set estimation mode to MCMC
EMODe 0 - set estimation mode to (R)IGLS
METH 1 - set estimation method to IGLS
EMODe 0 - set estimation mode to (R)IGLS
METH 1 - set estimation method to IGLS

[model specification]

START - begin (R)IGLS estimation

2. For MCMC:

tole 2 - set (R)IGLS convergence tolerance to 2
MTOT 5000 - set desired length of MCMC chain to 5000
EMODe 3 - set estimation mode to MCMC
EMODe 3 - set estimation mode to MCMC

BATCh 1 - allow start/next command to run for more than one iteration
MAXI - set maximum number of (R)IGLS iterations - note this is missing the value
START - begin (R)IGLS estimation

Below is an example of setting up and running a simple model using the tutorial data provided with MLwiN:

model setup
RESP "normexam" - set response to normexam
IDEN 2 "school" - set level-2 identifier to school
IDEN 1 "student" - set level-1 identifier to student
ADDT "cons" - add constant term
SETV 2 "cons" - set constant term to be random at level-2
SETV 1 "cons" - set constant term to be random at level-1
ADDT "standlrt" - add standlrt term

IGLS estimation
METH 1 - set estimation method of IGLS (this is the default for a new worksheet)
MAXI 20 - allow up to 20 iterations (this is the default for a new worksheet)
BATC 1 - allow more than one iteration when using start/next commands
STAR - begin estimation

MCMC estimation (note that unlike using the GUI I am not using residual starting values in order to simplify the commands)
EMOD 3 - let the GUI know that we are about to fit a model with MCMC
MCMC 0 500 1 5.8 50 10 1 1 1 1 1 1 - MCMC burnin for 500 iterations
MCMC 1 5000 1 c1090 c1091 c1003 c1004 1 1 - MCMC iteration for 5000 iterations
PUPN c1003 c1004 - summarise MCMC chains into C1096-C1099



If the response type for your model is different from the above, or you want to use different MCMC settings then you will need to modify the MCMC commands above to reflect this. If you have access to R or Stata then I would suggest looking at the scripts that are generated by R2MLwiN (see http://www.bristol.ac.uk/cmm/software/r2mlwin/) or runmlwin (see http://www.bristol.ac.uk/cmm/software/runmlwin/) as this will give you a quick way to see how to adapt the commands for different model types.
yongjookim78
Posts: 45
Joined: Tue Jan 10, 2017 3:36 am

Re: MCMC Macro

Post by yongjookim78 »

Dear Chris,

This is extremely helpful and I really appreciate it!!

Your codes worked for both (four level) linear and logistic models.

For instance, for linear model, I used this approach for model set up:

RESP 'linear_outcome'
IDEN 4 'region'
IDEN 3 'psu'
IDEN 2 'household'
IDEN 1 'id'
RDISt 1 3
RCON
ADDT 'cons'
SETV 4 'cons'
SETV 3 'cons'
SETV 2 'cons'
SETV 1 'cons'
ESTM 1
WSET
EXPA 2
WSET
ESTM 2
WSET
EXPA 2
WSET

For logistic model,

RESP 'binary_outcome'
IDEN 4 'region'
IDEN 3 'psu'
IDEN 2 'household
IDEN 1 'id'
RDISt 1 0
LFUN 0
DOFFs 1 'cons'
ADDT 'cons'
SETV 4 'cons'
SETV 3 'cons'
SETV 2 'cons'
SETV 1 'cons'
ESTM 1
WSET
EXPA 2
WSET
ESTM 2
WSET
EXPA 2
WSET

(For model set up, I put whatever the command interface showed when I did the same by using manual (mouse click) approach.)

Then, for IGLS estimation:

METH 1
MAXI 20
BATC 1
START

And for MCMC estimation:

EMOD 3
MCMC 0 500 1 5.8 50 10 1 1 1 1 1 1
MCMC 1 5000 1 c1090 c1091 c1003 c1004 1 1
PUPN c1003 c1004

My follow-up procedures are to save each model specification and compare them (i.e., point estimate, CI, p-value):

MSTO '\Null_Model'
STOR

I found that for both linear and logistic models, your codes worked well!

Here, I have got a few more questions.

First, for MCMC GUI, can I keep using the MCMC code that you provided once I make appropriate modification in the model set up part for different types of response variable?

Second, may I also need to modify "c1090 c1091 c1003 c1004" and "c1003 c1004" if I want to save and compare models?

Thank you very much again!

Best,
YongJoo
ChrisCharlton
Posts: 1250
Joined: Mon Oct 19, 2009 10:34 am

Re: MCMC Macro

Post by ChrisCharlton »

In theory you should adjust the MCMC 0 command to use metropolis-hastings when fitting a logistic model, however in practice the command will notice that you have chosen an invalid option and will use a default method. You will however need to change the specified model type (the last number in the MCMC 0 and MCMC 1 commands from one to two, to indicate that the model should be treated as binomial instead of normal. The command for the logit model should therefore be something like (I have put the changed values in bold):

MCMC 0 500 1 5.8 50 10 2 2 2 1 1 2
MCMC 1 5000 1 c1090 c1091 c1003 c1004 1 2

The columns that you mention are used as follows:
C1090 - chain of parameter estimates
C1091 - chain of deviance values
C1003 - parameter estimates
C1004 - parameter variances

Unless you need these to manually calculate additional statistics that are not included in the model comparison window then you can keep reusing the same columns between models, and these will be cleared when the new model is estimated.
yongjookim78
Posts: 45
Joined: Tue Jan 10, 2017 3:36 am

Re: MCMC Macro

Post by yongjookim78 »

Dear Chris,

Again, this is really useful and helpful, and I appreciate it a lot!

By the way, I have encountered another issue.

Until when I built a four-level null linear model and stored it, all worked well!

Then, when I added a set of variables, ran IGLS and MCMC, and stored it using a different model name, my MLwiN 2.36 as well as MLwiN 3.00 beta stopped working and the program was just closed.

My second codes are as below, and would appreciate if you could help me again.

After doing the prior procedure by using the codes that you kindly provided:

- Model set up:
EMODe 3
EMODe 0
METH 1
EMODe 0
METH 1

CENT 0
ADDT 'age'
CENT 0
ADDT 'female' 0
CENT 0
ADDT 'edu_2cat' 2
CENT 0
ADDT 'marri_3cat' 1
CENT 0
ADDT 'ho_income_3cat' 3
CENT 0
ADDT 'urbanicity' 1

- IGLS:
METH 1
MAXI 20
BATC 1
START

- MCMC:
EMOD 3
MCMC 0 500 1 5.8 50 10 1 1 1 1 1 1
MCMC 1 5000 1 c1090 c1091 c1003 c1004 1 1
PUPN c1003 c1004

- storing model and saving the work:
MSTO '\M1'
STOR

Thank you very much again!!

Best,
YongJoo
ChrisCharlton
Posts: 1250
Joined: Mon Oct 19, 2009 10:34 am

Re: MCMC Macro

Post by ChrisCharlton »

I can't see anything obviously wrong with the syntax. If you are able to replicate the problem using one of the MLwiN example datasets than I can run it in a debugger to see exactly where it is going wrong. Otherwise I suggest running each line via the command interface window to determine which command is causing the crash.
yongjookim78
Posts: 45
Joined: Tue Jan 10, 2017 3:36 am

Re: MCMC Macro

Post by yongjookim78 »

Hi Chris,

I just found that the codes actually worked when I ran the models having the data files on C: drive. The issue happened when I ran the models having the data files on flash drive. So, now it worked well, and I really appreciate all your great advices!

By the way, I have just got another quick question. For the MCMC codes, you kindly let me know that:

- For linear model:
MCMC 0 500 1 5.8 50 10 1 1 1 1 1 1
MCMC 1 5000 1 c1090 c1091 c1003 c1004 1 1

- For logistic model:
MCMC 0 500 1 5.8 50 10 2 2 2 1 1 2
MCMC 1 5000 1 c1090 c1091 c1003 c1004 1 2

Then, could you teach me how the codes need to be modified for 1) poisson model (i.e., to estimate binary outcome with prevalence around 20% with log link function), and 2) multinomial model (i.e., to estimate categorical outcome with three response options).

Thank you so much again!!

Best wishes,
YongJoo
ChrisCharlton
Posts: 1250
Joined: Mon Oct 19, 2009 10:34 am

Re: MCMC Macro

Post by ChrisCharlton »

As is partially documented in the help file for the MCMC command the codes for different distribution types (i.e. the last value in the MCMC command) are as follows:

1 - Normal
2 - Binomial
3 - Poisson
4 - Multivariate Normal
5 - Multivariate Mixed
6 - Unordered Multinomial
7 - Ordered Multinomial
8 - Negative Binomial
Post Reply