Recently I need to work on a very complicated hierarchical nonlinear model, the estimation requires complex estimation procedures. Therefore I resorted to the so called MCMC technique (Markov Chain Monte Carlo).
Winbugs is very powerful, and I absolutely like it. The problem with winbugs is that I need to use several softwares and go through a complicated process to test different models, different initial values, and the parameters to be estimated goes from 160 to more than 300.
First I create my regression datasets in MySQL (sometimes with a little help in ACCESS), then use a PERL script to turn them into winbugs format (I can post this program some time later, it is too ugly at this stage). After that I click the buttons in winbugs to do the study, if you ever used winbugs, you know the pain of clicking them. Due to the large number of parameters, the iterations are quite time consuming. I’m a cautious person, so I usually run 10,000 burn-ins, then use 50,000 iterations to calculate the estimates. After that, I have to select the text in the log window and paste the result to Excel, clean up the format and import to ACCESS or MySQL. Match the result with my regression covariates and export to SAS. In SAS, I do some simple regression and bootstrap analysis of the result. Usually, I will have some idea to test the next model, then I have to go over the whole process again from the beginning.
If you think this is a pain, you should see how I first used winbugs, I did not know how to estimate the 160 parameters, so I run 80 times of the model, each time obtain two values.
Today, I finally decided to put an end to this manual way of doing research. So I searched ways to do it combining WINBUGS and SAS. There are quite a few useful links on WINBUGS website regarding SAS. But they are relatively independent. I hope to write this small article to provide a unified way of using WINBUGS with SAS.
First of all, to automate the estimation, I need to be able to convert SAS datasets to winbugs format then submit my models within SAS, and without having to click on the winbugs window. The results should be read back to SAS.
The best candidate to submit the winbugs command, and gets it translated to a script is provided here by S.M. Mwalili. The documentation is given here. The idea is to submit a command of the following format:
%SasBugs(Work = c:/BugsSAS/schools/, Chains = 3, Burn_in =3000, Updates = 9000, Parameters = mu,sigma.theta theta, Coda = T, Dens = F, debug = F);
The problem with this solution is that you have to have everything (the model, the data, the init values) readily saved in winbugs format in that WORK directory. The model can be put there easily, but it still is a headache to convert and generate the data file. Another issue is that this procedure does not help to read back the results.
Rodney Sparapani’s and Matt Hayat’s SAS macros for BUGS
To convert between SAS and winbugs data format, we need some other macros available here. (Documentation)
_LEXPORT: Export a SAS dataset to a “list” data file
CODA2SAS: Converting CODA files to a SAS dataset
BUGStoSAS: converting BUGS output files to a sas dataset (Documentation)
SAStoCODA: Converting a sas dataset to CODA input files (Documentation)
- unzip them into a directory that I’ll call SASMACRO, but you can call it anything that you want
- add this line to your SAS configuration file, sasv8.cfg or sasv9.cfg
-insert sasautos ‘SASMACRO’
- several people have written to me about their broken SAS installations; if SAS can’t find it’s own SAS macros, like LOWCASE, try adding this (Unix):
-insert sasautos ‘!SASROOT/sasautos’
or this (Windows):
-insert sasautos ‘!SASROOT\core\sasmacro’
- now, call the SAS macros from your SAS program like:
%_decoda(out=dataset, chains=2, infile=coda);
%_limport(infile=alc.tbl, file=import.sas, out=alc);
%_lexport(file=alc.txt, data=alc, var=status alcgrams);
%_sexport(data=matrix, var=col1-col5, file=m.struct);
Zhiyong Zhang’s SAS WINBUGS pages
Zhiyong Zhang at Virginia Tech improved the Rodney’s macros and made them available here.
Two files got modified. The first one is _sexport.sas. It was changed to output the sample size (N) and number of variables (T) in the WinBUGS format data file. The second one is coda2sas.sas. It was modified to generate the history plot for each parameter.
Here gives a few links at to other resources here.
A good summary about how to implement the whole thing is in a PDF file here.