@Program name: block.prg @ 
@Purpose: to conduct bootstrap block or iid re-sampling or cases or residuals @ 
@Used for: 1996 course in non-parametric inference, Essex Summer School @ 
@Author: Chris Mooney @ 
@Date: June 20, 1996 @ 

new; 
t1=hsec;@for setting the timer@ 
output file = reg_res.out; @Defining output file@ 
output reset; @Overwrite previous output file @ 
library pgraph; 
graphset; 
@rndseed 14;@ 

@Set parameters @ 
resamp = 1; @set regression re-sampling method: 1 = residuals, 2 = cases@ 
re_block = 2; @set dependent data re-sampling method 1 = iid, 2 = block@ 
n=30; @sample size @ 
b=1000; @number of bootstrap re-samples @ 
r_weight = 5; @set random error weight @ 
s_weight = 30; @set sequential error weight@ 
t_cal=2.101; @set t-value for OLS CI. NOTE:needs to be reset with sample size@ 
block=5 ; @set block size@ 
k=n/block; @k is the number of blocks to be re-sampled@ 
tbeta = 2.0; tcons = 1.0; @pseudo-pop. parameters for regression data@ 

bb2 = zeros(B,1); @set up empty vector to hold bootstrapped beta-hat's@ 
i=1; 

@Generating pseudo-population data @ 
x1=ones(n,1); @generating vector of ones for constant term @ 
@x2=seqa(1,1,n);@ @generating fixed independent variable @ 
x2=rndu(n,1) * n; @generating random independent variable@ 
x=x1~x2; @set up X data matrix@ 
e = seqa(-1,2/n,n) * s_weight; @Generating sequential error @ 
e_r = r_weight * rndn(n,1); @Generating random error component@ 
y=(tbeta*x2) + tcons+ (e+e_r) ; @Generating the sample dependent variable@ 

__output=0 ; @to reduce superfluous output @ 
_olsres=1; @orders residuals to be saved @ 
{vn,m,bet,sb,vc,sd,si,cx,r2,resid,dw}=ols(0,y,x);@Regress full sample Y on X@ 
resid= y - (bet[1,.] + bet[2,.]*x2); @ save residuals for resamples@ 

print; 
print "sample size =" n; 
print "number of bootstrap re-smples = " b; 
print "regression re-sampling method (1=residual, 2=cases)=" resamp; 
print "dependent data re-sampling method (1=iid, 2=block)=" re_block; 
print "block size = " block; 
print; 
print "Full sample and OLS results:"; 
print "R-squared =" r2; 
print "Durbin-Watson statistic =" dw; 
print "OLS beta =" bet[2,.]; 
print "OLS se(beta) =" sd[2,.]; 
print "OLS confidence interval = (" (bet[2,.] - (sd[2,.]*t_cal)) "," (bet[2,.] 
+ (sd[2,.]*t_cal)) ")"; 
print; 

@ RESAMPLING@ 
do while i<=b; @Starting the re-sampling loop@ 
if re_block ==1; @Responding to iid resampling flag@ 
cindex=ceil(rndu(N,1)*N); @Setting up cases re-sampling index@ 
xc=x[cindex,.];yc=y[cindex,1]; @Re-sampling cases @ 
re=submat(resid,ceil(rndu(N,1)*N)',0); @Re-sampling residuals@ 
yr=x*bet+re; @Set up re-sampled Y vector for residuals re-sampling@ 

else; @Responding to block re-sampling flag@ 
nj=1; 
do while nj<=k; @Begin the block index development loop@ 
rn=(ceil(rndu(1,1)*n)); @rn is a randomly selected 
integer between 1 and n, to 
be used as the starting 
point for an index block@ 
t=seqa(rn,1,block); @Defines an index block starting from 
rn and increasing sequentially 
by 1's, of block length@ 
if nj==1; @Set t2 equal to the first index block@ 
t2=t; 
else; 
t2=t2|t; @Concatenate subsequent index blocks onto t2 vertically@ 
endif; 

nj=nj+1; @Increase block re-sampling loop index@ 
endo; @end the index development loop@ 

xcon = x|x; @Concatenating x, y, and resids on themselves@ 
ycon = y|y; 
residcon = resid|resid; 
xc=xcon[t2,.];yc=ycon[t2,1]; @Re-sampling cases @ 
re=residcon[t2,1]; @Re-sampling residuals@ 
yr=x*bet+re; @Set up re-sampled Y vector for residuals re-sampling@ 

endif; @end of block re-sampling flag response@ 

if resamp==1; @Responding to re-sampling flag@ 
{b11,se,r2}= regress(x,yr); @Residuals re-samplng regression@ 
else; @Responding to re-sampling flag@ 
{b11,se,r2}=regress(xc,yc); @Cases re-sampling regression@ 
endif; @end of loop choosing which type of re-sampling@ 

bsb2=b11[2,.]; @Identify the regression coefficient@ 
bb2[i,1] = bsb2; @Set up vector of bootstrapped beta-hat*'s@ 

i=i+1; @Adds 1 to bootstrap looping index@ 
endo; @end of resampling loop@ 

boot_b=meanc(bb2); @Gives the average of the beta-hat*'s, the 
bootstrap estimate of @ 
sd_b=stdc(bb2); @Gives the standard deviation of the beta-hat*'s, 
the bootstrap estimate of @ 
bb2=sorthc(bb2,1); @Sorting the vector of beta-hat*'s @ 
ll = ceil((b*0.05)/2); @Setting up lower and upper 
percentile points at .025 of b @ 
ul= ceil((b-ll)+1); 
low_b = bb2[ll,.]; @Gives the value of beta-hat* at the lower 
percentile point set by ll @ 
up_b = bb2[ul,.]; @Gives the value of beta-hat* at the upper 
percentile point set by ul @ 

print "Block re-sampling bootstrap results:"; 
print "Mean b* = " boot_b; 
print "Sd of b* =" sd_b; 
print "2.5th percentile point of b* = " low_b; 
print "97.5th percentile point of b* = " up_b; 
print; print "Minutes elapsed = " (hsec-t1)/6000; 
print; 

@REGRESS proc included to perform OLS@ 
proc(3) = regress (x,y); 
local xxi,beta,e,sse,sd,t,r2,A,il,n1; 
n1=rows(x); 
XXi=inv(X'X); 
beta=XXi*(X'Y); 
e=Y-X*beta; 
sse=e'e/(rows(X)-cols(x)); 
sd=sqrt(diag(sse*XXi)); 
t=beta./sd ; 
il=ones(n1,1); 
a=(eye(N1)-il*il'/n1); 
r2= (1-e'e/(y'*A*y)); 
retp(beta,sd,R2); 
output off; 
endp;