www.cmm.bristol.ac.uk/forum

I am trying to estimate a cross-classified random intercepts ordered multinomial model with common coefficients for predictor variables. Although I found no examples of such a model, I am assuming that MLwiN can run such models. (I know that MLwiN can run a cross-classified random intercepts binomial model.)

Yet I am having some problems sending the correctly specified model to MLwinN.

Starting with "xc.dta" (http://www.bristol.ac.uk/cmm/media/runmlwin/xc, clear) from "help runmwlin", I first modify "pupil" so that it is unique (apparently required for multimonial models in MLwinN).

.gen pupil2 = _n

I then generate the variables needed for the cross-classification: the pseudo-level 3 variables "s1", "s2", . . . "s19" to replace variable "sid".

.tabulate sid, generate(s)

I then generate the constraints that set the covariances of the pseudo-level 3 variables "s1", "s2", . . . "s19" equal to zero. (Note that I am using attain as an ordinal variable - there are 10 levels and, therefore, 9 contrasts.)

.local c = 1
.forvalues s=2/19 {
constraint define `c' [RP3]var(s1_123456789) = [RP3]var(s`s'_123456789)
local c = `c' + 1
}

When I run the following command . . .

. runmlwin attain cons, ///
level3(cons: (s1-s19, contrast(1/9)), diagonal) ///
level2(pid: (cons, contrast(1/9))) ///
level1(pupil2) ///
discrete( ///
distribution(multinomial) ///
link(ologit) ///
denominator(cons) ///
basecategory(10) ///
) ///
constraints(1/18)

I get the following error "detected my MLN":

"error while obeying batch file C:\~Temp\ST_0x00007.tmp at line number 100: FPAR 0 's1'
undefined explanatory variable referenced"

I suspect that my syntax for the (pseudo) level 3 is incorrect, but I can't figure out why. Thank you for any advice or help.

BP Lucas

Hi BP,

There are a couple of important things to note straight away

(1) MLwiN can fit discrete response models either by IGLS using quasilikeihood methods (MQL1, MQL2, PQL1, and PQL2) or by MCMC methods. Quasilikelihood methods are fast but approximate (i.e. biased) and the bias can be appreciable in certain situations (when you have small clusters and when the overall response proportion is extreme). Quasilikelihood methods are therefore not suitable for any final model or for publication purposes. MCMC does not suffer from these problems and the MCMC estimates will effectively be equal to unbiased maximum likelihood estimates.

(2) MLwiN can fit cross-classified models either by IGLS or by MCMC. The IGLS algorithm for fitting cross-classifed models are computationaly very inefficient (i.e. cumbersome and slow, only work for small data sets). The MCMC algorithm is computationaly much more efficient (i.e. much faster and works for large data sets and models with many classifications).

The above are two very strong reasons for using MCMC methods for fitting your ordianl response cross-classified model.

I suggest you read chapters 13 and 15 of the MLwiN MCMC Manual and then read and run the runmlwin do-files to replicate all the analyses presented in those chapters

http://www.bristol.ac.uk/cmm/media/runm ... sponses.do
http://www.bristol.ac.uk/cmm/media/runm ... _Models.do

And then have a go at setting up your model.

Start simple, first set up a binary response cross-classified model. Then first a three category ordinal response model and so on.

Best wishes

George

George - Thanks for the quick response.

I should have mentioned my reason for using likelihood-based methods. There are 2 reasons.

First, I want better starting values for MCMC. When I attempt to fit my models using MCMC, I have computational problems. For example, I get the following message "Warning: variance matrix is nonsymmetric or highly singular." I have tried to fix this by increasing the burnin length (up to 50,000) and tweaking the starting values.

Draper suggests using likelihood-based estimates as starting values for MCMC methods (Draper D. Bayseian multilevel analysis and MCMC. In: de Leeuw J, Meijer E, eds. Handbook of multilevel analysis. Springer, New York, NY; 2008). When likelihood-based estimates are not available, he suggests setting fixed-effects coefficients to zero and random-effects coefficients to one. Yet he mentions that "this approach may fail in the sense that the MCMC software is unable to begin sampling from such a poor starting place" when the posterior is "multimodal".

Perhaps my model cannot be fit from my data. Before I conclude that, however, I wanted to use the best starting values, which I was hoping to get from likelihood-based methods.

Second, in case my (theoretical) model can not be fitted with my data, I wanted to use likelihood-based methods to build an empiric model. As Goldstein says, "For purposes of data driven model selection, MCMC methods can be very time consuming. Generally therefore, we can use likelihood based methods to choose one or two suitable models that fit the data and then carry out the MCMC estimation" (Goldstein H. Multilevel statistical models. Wiley & Sons, West Sussex, UK;2011, p 54.)

Perhaps in my case (cross-classified ordinal model), however, it will still make more sense to use MCMC methods, because as you say, the likelihood-based methods are cumbersome and inefficient (as I have found).

Incidentally, most of the code from my original post was from the runmlwin do-files that try to replicate the findings in chapters 13 and 15 of MCMC manual. I think that I understand and can implement the material in each chapter. My difficulty is more with trying to *combine* the material from both chapters. If you have any runmlwin syntax from a cross-classified ordinal model, I'd love to see it - if for nothing else, just to see where I have been getting it wrong.

Otherwise, any additional advice would be greatly appreciated.

Thanks.

BPL

Hi BPL,

Typically the fixed part parameter estimates are not that sensitive to the random part model specification. You can therefore do most of your fixed part model exploration by simply fitting two-level ordinal models using IGLS. Once you have a couple of suitable final models you can then move to MCMC to introduce the cross-classified random part structutre. This also means that fitting the two-level ordinal model (by IGLS and then MCMC) should give you fairly decent starting values for the fixed part parameters of your cross-classified ordinal model.

You can retrieve these parameter estimates and store them as a vector by typing
As your next moel will be the cross-classified ordinal model, it will have one extra parameter (a new variance associated with the new classification).

You will therefore need to make this vector longer by one element. Be careful here, you will need to insert the new parameter in the correct place within the vector a.

For a sensible starting value, I suggest you just set each of the two variance parameters in your cross-classified model to be half the estimate of the single variance parameter in your two-level ordinal model. Leave all fixed part parameters estimates as they were in the two-level ordinal model.

In MLwiN you can get the message "Warning: variance matrix is nonsymmetric or highly singular." when one of the random part parameter matrices has negative variances or covariance parameters which imply correlations which lie outside (-1,1). IGLS can give these nonsensical estimates when the model is poorly specified or has had convergence difficulties. To fit the model by MCMC you need to manually edit the starting values so that they are well behaved (i.e. positive variances and correlations between (-1,1)).

Best wishes

George