***
*** Example 2.1. Comparing OLS, Tobit, and Powell's censored LAD estimators
*** 
*** The heteroskedastic censored regression model is discussed on 
*** pp. 85-90. I give the code for the Monte Carlo experiment 
*** comparing Powell's censored LAD estimator with OLS and MLE 
*** Tobit, the results of which are discussed on p. 90.
*** 
*** Note: The following provides the STATA code used to produce 
*** results used in "The Analysis of Household Surveys: 
*** A Microeconometric Approach to Development Policy", by Angus Deaton.
*** This book, published for the World Bank by The Johns Hopkins University
*** Press and scheduled for release in 1997.
***
version 4.0
#delimit ;
set more 1; 

/* control program */ 

/* number of observations */ 

/* either 100 or 1000 */
global nobs=100; 

cap program drop doit;
program define doit;
local ct = 1;
while `ct' < `1' {;
dis ("Experiment number `ct'");
runtest;
append using temp;
qui save temp, replace;
local ct=`ct'+1;
};
end; 

cap program drop runtest;
program define runtest;
drop _all;
qui set obs $nobs;
gen x=_n/($nobs/100);
gen u=invnorm(uniform());
qui replace u=u*20;
gen ytrue=max(0,x-40);
/* this is the heteroskedasticity */
qui replace u=u*(1+0.03*ytrue);
gen ystar=-40+x+u;
qui replace ystar=0 if ystar < 0;
qui regress ystar x;
/* ols regression */
global bols=_coef[x];
qui tobit ystar x, ll;
/* mle tobit estimation */
mac def btob=_coef[x];
/* censored LAD estimation: 10 loops */
qui qreg ystar x, qu(0.5);
predict xb1;
local iloop=2;
while `iloop' <=10 {;
qui qreg ystar x, qu(0.5), if xb1 > 0;
predict xb2;
qui replace xb1=xb2;
drop xb2;
local iloop=`iloop'+1;
};
/* powell estimate coefficient */
mac def bpow=_coef[x];
drop _all;
qui set obs 1;
gen zols=$bols;
gen ztob=$btob;
gen zpow=$bpow;
list;
end; 

drop _all;
set obs 1;
gen zols=999;
gen ztob=999;
gen zpow=999;
cap save temp, replace;
doit 100;
drop if zols==999;
gen absols=abs(zols-1);
gen abstob=abs(ztob-1);
gen abspow=abs(zpow-1);
gen powbest=0;
replace powbest=1 if (abspow <=abstob) & (abspow < absols); 

summ; 

[Return to Program List]