MCMC Error message?

[GDistiller]

Hi
I am trying to run a multistate model with random effects. I expect it to take at least a few days to run (started it yesterday) but when I came in this morning I saw a whole lot of error messages in the Mark run window. These error messages say:

"* * Warning * * Error number 43 returned from DPOCO in MCMC"

This seems to apply to certain iterations so for example even though I have a whole screen full of these messages I can see that it is not happening for every iteration. For example it says: Model name Series 1, Iteration 130 followed by the error message and then on the next line it says: Model name Series 1, Iteration 139 followed by the error message.

I'm not sure how serious this is and if I should just let it run? Or if I should abort the process? Any advice as to what error 43 means and how it can be resolved would also be greatly appreciated.

Thanks!

Greg

[gwhite]

Greg:
I'm guessing that you have specified a VC matrix, and that this matrix is not always positive-definite. You probably need to re-think what all the covariances have to look like given the ones you specified. If you didn't specify a VC matrix, then I'm not sure how you are getting this error.
Gary

[GDistiller]

Thanks Gary...

I did specify a VC matrix with correlations between my different hyperdistributions for survival. I'm not sure what you mean by "You probably need to re-think what all the covariances have to look like given the ones you specified." I have 3 physical states (locations) and want to estimate the correlation in survival between the states. I have followed the example in the help so that I correlate survival at time i in state 1 with survival at time i in state 2 etc. In total I am therefore estimating 3 rhos, survival between state 1 and 2, 1 and 3, and 2 and 3...

Do you think I must have made a mistake in the VC matrix?

Thanks very much!

Greg

[gwhite]

Greg:
I suspect you did make a mistake, because as you describe he matrix, it should be positive-definite.

Gary

[GDistiller]

Hi Gary
OK thanks, I thought I had set it up carefully but will relook at what I did and try get it to run.

Do you mind if I ask you another question? Its related to a seperate posting I made that I have not had any replies to, I've been thinking about this quite a bit and not managed to get any advice yet so would really appreciate any help. The question concerms how to parameterize the MCMC model when one has dummy coded indicator variables for the time / years, and specifically I am talking about the intercept. My confusion stems from the fact that the intercept term is actually giving the estimate for the reference year, yet the few examples I have found seem to exclude this from the temporal hyperdistribution. My understanding of this is that the intercept would then estimate the mean survival and the random effects would be random deviations around this mean (and so one would set the mean of the hyperdistribution to be zero). I'm not sure if this is correct as conceptually it seems like all yearly estimates (incl the reference year) should come from the same dbn? But at the same time it also seems odd to me to not have an intercept as part of the fixed part of the model.

I then started thinking about a std linear mixed model where one has y=XB+Zu. Typically the intercept would be part of the fixed part and the Z design matrix would include a random effect for every year ie if there are t years then there would be t columns in Z so that it differs from the fixed effect coding where you have t-1 columns to code the t years. So I started thinking that I should setup a design matrix in Mark that has an intercept as well as t indicator variables for my t years. I would then include these t parameters into a common hyperdistribution where again the mean is set to be zero. This way the intercept would be an overall mean and every year in the data would be able to have a random effect included. One problem with this though is that I would imagine that I cannot run it first (before the MCMC run) due to the redundancy in parameters. I could of course set it up and provide starting values without running it as a fixed effects model first. Does this sound correct or is it better to just exclude the reference year from the hyperdistribution?

Many thanks for the help...

Greg

[gwhite]

Greg:
I'm assuming you only have a time effect, so code t without an intercept. Then, you can specify a hyperdistribution across all of the t effects.

If you require an intercept (e.g., where you have both group and t effects), then you can specify the t-1 effects using the -1 trick for the last effect. Then, the intercept is the mean of the t effects, not the reference effect. So

1 1 0 0 0
1 0 1 0 0
1 0 0 1 0
1 0 0 0 1
1 -1 -1 -1 -1

you specify the hyperdistribution on columns 2-5, and force the mean to be zero via the prior.

Gary

[GDistiller]

Hi Gary
OK thanks very much for the reply. I do have a group effect as well (the effect of being a transient) so will follow your advice on how to set it up. Would the intercept then be the mean of the t effects but in a particular group (say for group=0)? I had thought of changing the dummy coding scheme but was not sure as I have come an example where the intercept is excluded from the hyperdbn but the dummy coding is just the regular ref category scheme (there is a group effect as well). Is that actually wrong?

Please can you clarify how I force the mean to be zero via the prior? I know that one gets the "Hyperdistribution parameters" dialog box and can supply 4 values for each hyperdistribution. If the mean of the prior (so the 3rd value) is zero does this effectively force the mean of the random effects to be zero?

Lastly, I have tried to rerun my model where I specify the VC matrix and was getting that MCMC error that is the topic of this initial posting. I have looked at the matrix I have specified carefully but just cannot see what I am doing wrong. I have also shown the matrix to a colleague and we have coded it in Mathematica to check if it is positive definite. It seems to be fine but for certain values it is not positive definite - but I thought you had said there must be a problem with the structure rather than actual estimates for certain iterations? Should I leave it to run and hope that after the burn in period the chain will have settled in a space where the estimates will result in the matrix being positive definite?

Many many thanks for all your help - I really do appreciate it!

Greg

[gwhite]

Greg:
You can use the -1 coding for the group effect as well so that the intercept is now the mean across all times and all groups on the logit scale.

To force the mean of the hyperdistribution to be zero (or very nearly so), set the mean of the prior distribution to zero and the sigma to 0.001 (or even try 0.0001), so that the prior is very narrow, and does not allow the mean to deviate from zero.

I'm not seeing how the VC matrix can not be positive definite if you are specifying the form that I think you are. However, as you are speculating, when the VC matrix is not positive definite, then that set of values is never accepted, and so is not included in the posterior distribution. If some range of a correlation is not possible, set the prior on the correlation to exclude the impossible values. The prior on a correlation defaults to a uniform [-1, 1]. Specify different endpoints.

Gary

[GDistiller]

Hi Gary
Got it re the group effect, should have realised that myself...

Wrt the VC matrix, when we coded it in Mathematica we were using covariances rather then correlations on the off-diagonals. For certain values the matrix was not positive definite ie when violating the Cauchy–Schwarz inequality (that the squared cov(X,Y)<=product of the variances). I'll try it again using a correlation matrix and specifying different correlation values in the valid range of -1 to 1.

I do not see how any particular valid correlation values are impossible for this structure. This may be a bit presumptious of me but I have drawn up the matrix I am specifying in excel - is there any chance you would be prepared to look at it and see if there are any obvious errors that I am missing? Is there a way to upload an attachment to the forum?

Thanks again

Greg

[cooch]

... I have drawn up the matrix I am specifying in excel - is there any chance you would be prepared to look at it and see if there are any obvious errors that I am missing? Is there a way to upload an attachment to the forum?

The Excel spreadsheet Greg mentions can be downloaded here.

As an aside - no, there is no way to upload an attachment to the forum - I've turned that feature 'off', since in my experience allowing uploads quickly gets out of control, and floods the server (and the generated emails when things are posted to the forum). In this day and age, most people have access to a server (there are *lots* of free hosting providers for small files out there) - if not, I can post things from my end if I think there may be general interest (as in this case).