Stata FAQ
How can I graphically compare OLS and BLUP results in Stata?

There are two ways you could go about this: 1) graph the shrinkage of OLS to BLUP (Best Linear Unbiased Predictor) and 2) graph both the fitted OLS regression lines and the BLUP fitted lines on the same graph.

Below is a do-file that does both of these for the imm23.dta that was used in the Kreft and de Leeuw Introduction to multilevel modeling.

The shrinkage graph has slopes on the vertical axis and intercepts on the horizontal axis. The tails of the arrows give the OLS values and the points give the BLUP values.

Note that there is less variability in the BLUP intercepts which also have flatter slopes than their OLS counterparts.

/* generate ols intercepts and slopes using statsby */

use http://www.ats.ucla.edu/stat/stata/examples/mlm_imm/imm23.dta, clear
statsby , by(schi): regress math homework
rename _b_homework slope
rename _b_cons intercept
label var slope "ols slope"
label var intercept "ols intercept"
sort schid
save t2, replace

/* generate blup intercepts and slopes using xtmixed */

use http://www.ats.ucla.edu/stat/stata/examples/mlm_imm/imm23.dta, clear
egen grp=group(schid)
quietly xtmixed math homework || schid: homework, variance covar(un)
predict u1 u0, reffects
generate s = _b[homework] + u1
generate i = _b[_cons] + u0
label var s "blup slope"
label var i "blup intercept"

by schid, sort: generate n = _n
merge schid using t2
list schid slope intercept s i if n==1, clean

/* graph one -- shrinkage graph */

twoway pcarrow slope intercept s i if n==1
more

/* graphs two & three -- graphs of regression lines */

generate h0 = 0
generate h7 = 7
generate o0 = intercept + slope*h0
generate o7 = intercept + slope*h7
generate b0 = i + s*h0
generate b7 = i + s*h7

/* graph is split into two parts due to the number of fitted lines being displayed */

twoway (pcspike  o0 h0 o7 h7 if n==1 & grp<12)  ///
       (pcspike b0 h0 b7 h7 if n==1 &  grp<12), ///
       legend(off) title("blue ols - red blup")
more
twoway (pcspike o0 h0 o7 h7 if n==1 &  grp>11)  ///
       (pcspike b0 h0 b7 h7 if n==1 &  grp>11), ///
       legend(off) title("blue ols - red blup")

Below is the output generated by the do-file.

       schid       slope   interc~t           s          i  
  1.    6053    2.091476   51.61236    2.167148    51.1254  
 45.    6327   -8.566667       63.8   -6.033969   61.11028  
 53.    6467    6.981482   41.24074    6.158613   41.90844  
 58.    7194   -2.795455   52.65152   -2.215304   51.89136  
 82.    7472   -3.553797   50.68354   -2.952869   50.28959  
105.    7474   -3.306123   59.77551   -2.600884   57.82726  
122.    7801   -3.069767   54.69556   -2.453294   53.71066  
144.    7829   -2.920123   49.01229   -2.777317   49.28519  
164.    7930    7.909091      38.75     7.39565   39.50748  
188.   24371           5      38.35    4.681487   39.36426  
208.   24725    5.592664   34.39382    5.351505   35.48435  
230.   25456   -4.718411   53.93863   -3.234086   52.85429  
252.   25642   -2.486056   49.25896   -1.397812   48.43803  
272.   26537    3.874317   46.20492    3.900892   45.72621  
288.   46417     4.40731   47.64752    4.299445     47.299  
311.   47583    3.376682   41.76413    3.198847   42.51505  
331.   54344    3.507732   35.22423     3.24509   36.74875  
350.   62821     1.09464   59.21022      1.3621   57.89467  
417.   68448     6.49631   36.05535    5.893169   37.65525  
438.   68493        5.86      38.52     5.31894     39.459  
459.   72080    5.094132   43.03178    4.821579   43.33112  
486.   72292    6.335052   37.71392     6.01733   38.44124  
506.   72991       5.575   43.77857    5.397511   43.62207