11-26-2015 08:53 PM - edited 11-26-2015 09:49 PM
I've tried a simple regression on both SAS and Winbugs, but SAS generated different estimates everytime I hit "run" even though I didn't change any values (just press "run" again and again), while Winbugs gave me quite stable results.
Here is my SAS code, would someone tell me anything I missed in the program? Or SAS MCMC itself is not functioning well?
Input ID $ Y X_mu X_sigma;
A 15 4 1.1
B 8 1.2 0.8
C 12 2.8 1
Proc MCMC Data=base nmc=10000 Autocorlag=10000 thin=1 nbi=10000 monitor = (intercept beta X);
Parms intercept beta X;
Model Y ~ normal(Intercept + beta * X , sd=1);
Prior Intercept ~ normal(0, sd=100);
Prior beta ~ normal(0, sd=100);
Prior X ~ lognormal(X_mu, sd= X_sigma);
For running five times, SAS generated estimates of beta as -1.3, 7.8, -2.3, 0.4, 2.3 (and kept changing...),
and Winbugs generated a stable range 0.70, 0.75, 0.70, 0.70, 0.70.
11-28-2015 06:21 AM
If you want the same results, use the SEED= option in the PROC MCMC statement:
proc MCMC Data=base SEED=12345 ...;
This sets the random number seed. By default, SEED=0, and PROC MCMC gets a random number seed from the time of day.
11-29-2015 09:56 PM - edited 11-29-2015 09:57 PM
Thank you, @Rick_SAS. The Seed option did work, and the results are the same now.
But I checked the posterior values for each iteration, found that the values were barely moved no matter the seed was set or not (the first 20 values are as belows.) I need to increase the thinning rate, say 10, or the nmc to 200,000 to maker sure the results are right, but then it takes a lot of running time. Is this normal for SAS?
|Iteration||intercept||beta||x||Log Prior Density||Log-Likelihood||Log Posterior|
11-30-2015 01:53 PM
You are not getting any mixing at all in the chain. Your profile plots must look very strange. It is probably related to the fact that you do not have a free residual variance parameter. You are fixing the standard deviation (variance) at a constant in the model. This usually is a parameter with a prior.
I don't know why you are treating X as a parameter rather than as a fixed predictor variable. I assume you have a good reason (perhaps, as a measurement-error model).
11-30-2015 04:21 PM
The problem is your treatment of X as a parameter. I tried this with some data, and MCMC cannot mix. If you just had a model with X_mu as the predictor variable, it runs fine. Of course, you may have a reason to treat the predictor as a variable, although I am not aware of a reason for doing so.
11-30-2015 09:03 PM
Thank you, @lvm. This helps a lot.
yeah, I was guessing the "X" was the reason. When I set X constant, everything looks fine. But since X is an estimation instead of an observation in my work, I need to assign a distribution for X. I may need to figure out a way to solve this problem.