libname data "C:\Documents and Settings\ejohnson\My Documents\MultilevelModels";

* Get the region specific means;

proc sort data=data.nmmaps_pm_coef;
by region;
run;

proc means data=data.nmmaps_pm_coef;
var beta;
by region;
run;

* Get the overall mean of the betas;
proc means data=data.nmmaps_pm_coef;
var beta;
run;

* Run an ANOVA to determine if heterogeneity across regions is greater than;
* within region heterogeneity;
proc anova data=data.nmmaps_pm_coef;
class region;
model beta = region;
means region / hovtest=bartletts;
run;

proc mixed data=data.nmmaps_pm_coef;
class region;
model beta = region;
random intercept / subject=region;
run;

* NOTE, this is the correct specification of the model,
but it does not converge!  Don't know whats happending here!
Sorry!;
proc nlmixed method=FIRO data=data.nmmaps_pm_coef;
model beta~normal(theta+bj,sigma2);
random bj ~ normal(0,tau2) subject=region;
predict bj out=ebayes;
run;